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SUMMARY AND INTRODUCTION 


During the past decade several excellent computer codes for analyses 
of complicated flow fields have been developed. NASA -Langley’ s SEAGULL 
program (Appendix A) that treats multiple shockwaves and contact surfaces, 
using a floating shock fitting technique, is one such code. Previously, calcula- 
tions performed with this code were limited to the ideal gas assumption. It is 
the purpose of this document to describe the extension of the code to handle 
real gas effects via the incorporation of a general finite rate chemistry and 
vibrational energy exchange package. 

This document and the associated computer code (modified SEAGULL) 
now provide NASA -Langley with a benchmark finite rate general chemistry 
nozzle and plume code. The modified code retains all of its original features 
plus the capability to treat chemical and vibrational nonequilibrium chemistry. 
The chemistry package is extremely general in nature, handling any chemical 
reaction or vibrational energy exchange mechanism as long as thermodynamic 
data and rate constants are available for all participating species. 
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1- SYMBOLS 


Symbol 


A 

b(z) 

c(z) 

C 

P 

C 

v 

f(s) 

f. 

1 

F. 

1 

G 

h 

H 

k 



P 

q 


r 

3t 

R 


Description 
defined as used 
lower duct wall ordinate 
upper duct wall ordinate 
specific heat at constant pressure 
specific heat at constant volume 
entropy production (Eq.4) 
weighting factors (Eq. 15) 
species mole/mass ratio 
defined by Eq. (50) 
static enthalpy 
total enthalpy 
reaction rate constant 
equilibrium constant 
natural log of pressure 
pressure 
total velocity 
radial coordinate 
universal gas constant 

universal gas constant divided by global 
molecular weight 



Symbol 


Description 


s 

entropy 

T 

temperature 

u 

radial velocity component 

w 

axial velocity component 

w. 

1 

net production rate of species i 

X 

transformed radial coordinate 

y 

axial coordinate 

z 

axial coordinate 

Greek 

a. 

i 

species mole fraction 

(3 

defined by Eq. (3) 

r 

ratio of specific heats 

v', V 
1 1 

stoichiometric coefficients 

p 

density 

+i 

species molecular weight 


global molecular weight 

Subscripts 

i 

refers to species i 

in. 

refers to initial conditions 


2 



2. DISCUSSION 


2. 1 BASELINE NOZZLE CODE - PROGRAM SEAGULL 

Program SEAGULL is a computer code for the numerical analysis of 
complex two-dimensional or axisymmetric supersonic inviscid flows of a 
perfect gas. The fundamental limitation of this program is that the com- 
ponent of the Mach number in the axial direction at any point in the flow field 
must remain supersonic. The program was primarily designed for the analy- 
sis of internal flows. It can compute the flow field produced by a single duct 
or several ducts that are merged. Although designed for internal flows, it 
can also compute jets and plumes. The program continuously monitors the 
flow field to detect the formation of shock waves. All discontinuities are 
treated explicitly, and all interactions are treated by a locally exact solution. 

A description of the method used is given in Appendix A. 

2.2 MODIFIED NOZZLE CODE 

The modified program Seagull is now a computer code for the numerical 
analysis of complex two-dimensional or axisymmetric supersonic inviscid 
flows with fully coupled nonequilibrium chemistry. The code retains all the 
salient features discussed in Section 2.1, plus the additional capability to treat 
frozen or finite rate chemically reacting flows. The generalized chemistry 
package included in the modified Seagull also permits the analysis of vibra- 
tional nonequilibrium energy exchanges. 

2.3 DEVELOPMENT OF EQUATIONS 

The basic equations governing the flow of gases with nonequilibrium 
reacting chemistry are the following: 
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V • pq = 0 Global Continuity 




V q * pq + Vp = 0 Momentum 


pq • VH = 0 Energy 


pq • VF^ - w. = 0 Species Continuity 


N 


a. 


P = p ffT ^ ^ Equation of State 

i=l 1 


(1 , 


H = h + ^V(q . q) Total Enthalpy 


w. 

i 


M 

±L (v!' . - v \ .) 

1 i = l 


N 


v\ 


N 


V 


k n f/ ,j - k n f" 1 >j 

J f=i f D ^=i “ 


General Species 
Production Rate 




where N species are involved in M reactions. 

These equations are expanded in a Cartesian coordinate system, scaled 
(i.e., nondimensionalized) with respect to some initial reference conditions, 
and are then solved in the modified SEAGULL code. 

To avoid unnecessarily complicating the discussion, the remainder of 
the derivation is shown for a two-dimensional case only. The code, however, 
can also handle axisymmetrics flows. 
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Scaled and expanded, the equations then become 




w 


(3 P + wu + A , wu P - aw - £(s) = 0 

^ z r I r r ^ 

P w z + -| u r + A i wuw r “ A 1 ^ uP r + A^ f(s) 


= 0 


u + — u + — — r P - 0 
z w r wy r 


w ^ 

wF. + uF. - 77 = 0 
i i P 

z r 


H + — H = 0 

z w r 




( 2 ) 


where the subscripts refer to differentiation with respect to the r, z coordinates, 


(3 , A and A^ are conveniently defined as follows: 

.2 T 


(3 = A.w 




y. - 1 

A - 1 - in - ■ 

A 1 “ 1 A. 


a, = r. c 

2 in p 




J 


(3) 


The function f(s) is the entropy production due to finite rate chemical re- 
actions . 

f(s) = (q. vs) chem (4) 

The general set of governing equations (Eq. (2)) are then recast in 
terms of the computational coordinates (X, Z) that are defined by 

X = (r -b(z))/(c(z) - b(z)) 

and (5) 

Z = z 

where b(z) and c(z) define the lower and upper walls of the duct or nozzle. 
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All the primary variables in the equation set are scaled with respect 

to the initial conditions. The pressure, density, temperature and molecular 

weight are non-dimens ionalized by p.^, P^ n > » iji^* respectively. Velocity is 

conveniently scaled by ^p and entropy is non-dimensionalized by the 

specific heat at constant volume C , i.e., S = ( [S] - [S.J )/C . 

^ 1 »■» ^ 


The specific heat at constant pressure, C^, is scaled with respect to 
its initial value, [Cpl= [Cp]/[C p . ]• Static enthalpy is non-dimensionalized by 


m 

0^ T^ n> and for consistency the total enthalpy equation then becomes 


m 


H T = h + 


[ft] . 2 , 2 , 

[C] (u + w ^ 


( 6 ) 


m 


where all values are non-dimensional except [!}fj and [c 



The finite rate chemistry parameters, of course, have dimensional char- 
acteristics. Interrelationships between the flow field and the finite rate re- 
actions are treated internally by temporarily dimensioning the variables ap- 
pearing in the species continuity equations. The reacting chemistry calculation 
is then carried out and the results are subsequently non-dimensionalized. 


2.4 SOLUTION TECHNIQUE 

A finite difference numerical scheme utilizing the explicit MacCormack 
operator is used to integrate the equations for conservation of global mass, 
momentum and energy. The species conservation equations are solved by 
an implicit technique similar to that developed by Moretti (Ref. 1). This 
method was chosen to provide good stability characteristics in the overall 
numerical scheme. 


The general solution technique is discussed adequately in Appendix A. 

The same methodology is retained in the modified SEAGULL with the finite rate 
reacting chemistry. The only exception is the species continuity equation set 
which is solved by the following technique. 
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A detailed description of the rate processes that occur in finite rate 
reacting flows requires that a myriad of mechanisms be considered to in- 
clude all the possible chemical and vibrational reactions of dissociation, 
formation, recombination, etc. 


All of these, however, can be treated with a very general formalism. 

In the form usually quoted in chemical kinetics (Ref. 2) the phenomenological 
law of mass action states that the rate of a reaction is proportional to the 
product of the concentrations of the reactants. Thus, for a general reaction 
of the form ^ 


E 

j=i 


V\ A 




j=l 


V\' A. 

J 3 


(?) 


the net rate of production w^ for any participating species for which the 
stoichiometric coefficients and v'S are not equal can then be written as 




N 

77 
1 = 1 


1 


N . 

™ 9 \ 

n f 

1 = 1 


( 8 ) 


Assuming small deviations from equilibrium, the forward and backward 
reaction rate constants, k f and k^, respectively, can be related to the con- 
centration equilibrium constant and to the pressure equilibrium constant as 
follow s: 

N 

k E O'! - V'{) 

= K c = K p ( ST) 14 (9) 

The significance of the pressure equilibrium constant is that it can be 
easily evaluated for any reaction using tabulated values of the equilibrium 
constant for formation from the elements. Values of are commonly tabu- 
lated in conjunction with specific heats, entropies and enthalpies as a function 
of temperature, and are available in general for most species. An equally 
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convenient method exists for determining K 
accompanying the reaction, i.e.. 


P 


from the change of free energy 


K = exp ( - AG /ST) (10) 

P 

where AG is the change in free energy during the reaction process. Free 
energy values are also available for most species in tabular form. This 
method is used to compute in Program SEAGULL. The JANNAF thermo- 
chemical tables (Ref. 4) are used as the source for obtaining input data for 
various chemical systems. AG is computed directly by taking the difference 
in free energy between the products and reactants. 


For reasons of computational speed and efficiency, the program con- 
tains explicit expressions, as obtained from Eq. (11), for the most commonly 
encountered reaction mechanisms. Twelve types of reaction mechanisms 
are considered as possible contributors to the calculation of the net rate of 
production, w^. 


Reaction Type 


(1,7) 

A 

+ 

B 


C + 

D 

(2,8) 

A 

+ 

B + M 


C + 

M 

(3,9) 

A 

+ 

B 

— * 

c + 

D + 

(4, 10) 

A 

+ 

B 


c 


(5, 11) 

A 

+ 

M 


c + 

D + 

(6, 12) 

A 

+ 

M 

~ * 

c + 

M 


Reaction types (7) through (12) correspond to reaction types (1) through 
(6), but proceed in the forward direction only. 
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To reduce roundoff and truncation errors, the forward and backward 
rates for each reaction are computed separately. All contributions to the 
molar rate of production of a given species are then computed and added 
algebraically to form matrix coefficients (discussed later). Since reaction 
types (7) through (14) proceed in the forward direction only, the second term 
on the right-hand side of Eq. (11) is disregarded in calculating the contribu- 
tions to the coefficient matrix. 

In reactions (2), (5) and (6) as well as in (8), (11) and (12), M denotes 
a third body species which can be specified. For these reactions the situa- 
tion often occurs where for various third bodies the respective rate constants 
differ only by a constant multiplier. These multipliers can be considered as 
third body efficiencies or weighting factors. If such a case is encountered, 
the third body species mole mass ratio F^ becomes effectively a fictitious 
mole mass ratio, consisting of the weighted sum over all those species 
having a nonzero weighting factor, i.e., 

F .. = V* f. F., (12 

M i M. 

i 1 

where F are the weighting factors. 


The forward rate constant, k^, is generally a function of temperature. 

It is most often expressed in Arrhenius form. Again, for speed and efficiency 
in computation, the rate constants are divided into five types: 


Rate Constant Type 


( 1 ) 

k f = A 
1 -N 

( 2 ) 

k f = AT 

( 3 ) 

k^ = A exp(B/ !RT) 

( 4 ) 

k f = AT" N exp(B/ 9RT) 

( 5 ) 

k f = AT“ N exp(B/ !HT M ) 


( 13 ) 
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The equations presented in this document provide a very general forma- 
lism for the evaluation of various rate processes. The specification of par- 
ticular systems and associated rate constant will be up to the program user. 

Consider now the general species continuity equation 

pq . VF t = w. (14) 

and making use of the foregoing discussion of the rate process we now pro- 
ceed to describe a calculational technique for determining the individual 
species composition on a point-by-point basis. The description of this 
process is substantially simplified if Eq. (11) is specialized to a particular 
reaction type, say number (7) from Eq. (14) which is a one-way, two-body 
reaction. 

A + B -**C + D (15) 


the net production rate for this process is 

w = -k f P 2 F A F B (16) 

and the species continuity equation for species B then becomes 

pq- vf b = -k f p 2 f A F b (17) 

This equation can readily be solved using finite difference techniques 
employing explicit relationships, such as Euler or more sophisticated schemes, 
such as Runge-Kutta. The step size for integrating this equation, however, is 
severely limited by stability criteria. It can be seen from Eq. (20) that the 
rate of change of a species along the streamline becomes increasingly larger 
as the flow speed is slowed, the density increased, or for fast reaction rates. 

In rocket engine problems, combinations of slow speeds, high densities and 

fast reaction rates (i.e., quasi-equilibrium) are quite common and integration 

- 8 

step sizes so small (i.e., < 10 meters) are encountered that the solution 
becomes impractical in terms of computation time. 
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For this reason, the technique described in Ref. 3 based on a lineariza- 
tion of the production rates was utilized. Writing Eq. (20) in finite difference 
form over an integration step from station n at z = z q to n+1 at z = z q + Az 


= F. 


k^Azp 


n+1 


n 


n+1 


8f b 

. u • n 

- — Az — * 

B n+1 W 9r 


(18) 


And evaluating all the species concentrations at the downstream point results 
in a set of simultaneous nonlinear algebraic equations. The lateral gradient 
8F B 

term — g is evaluated at the upstream station n and uses windward differences. 

n 

In order to solve these equations we must then linearize the term F. F„ 

n+1 n+ 

which is accomplished following the lead of Ref. 3. If this term is expanded 
in terms of its values at station n along with the increments over n to n+1 
we can obtain the following expression. 


n+1 


B 


= F 


n+1 


n 


B 


+ F. 


n+1 


n 


n+1 



(19) 


neglecting products of differentials which are assumed to be of second order 
importance. Equation (21) can now be written in its linearized form. Let 
C = Az k^p/w, then 



and similarly. 



k f b a1 + f b f a „ - f a f b 1 

L n n+1 n n+1 n n-l 

[ f a f b +f b f a - f a f b 1 

L n+ 1 n n+ In n n J 


3F. 

u A ] 

— Az -z — 
w dr 


n 


u 

w 


Az 


3F 

FF 


n 


( 20 ) 
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Equation (20) can then be expressed in terms of a set of unknowns and calcu- 
lable coefficients, C. Rewriting these we obtain 


B + , = Q B - CF A < F B > - CF B < F A > 

n+1 n n n+1 n n+1 

'a +1 = q a - cf b < f a > - cf a < f b > 

n+1 n n n+1 n n+1 


( 21 ) 


where 


9F. 

Q. = F. +CF. F. - - Az 
i i i j w Or 

n n n J n n 


( 22 ) 


A (1 + c F b ) + ( c F a ) F b - Q a 

n+1 n n n+1 n 

’ (C F ) + (1 + C F ) F = Q 

A n+1 B n A n B n+1 B n 


(23) 


A matrix can now be formed using totally known information. 


1 + C F 


B 


C F 


n 


n 


C F, 


1 + C F 


n 


n 


n+1 


n+1 


Q 


Q, 


n 


n 


(24) 


The matrix [a] [x] = [bJ is then solved for the unknown compositions F A , 

A , 

F -,3 via a triangulation technique. Although consuming more time per 
B n+1. 

integration step than an explicit formulation, the implicit technique employed 
here is unconditionally stable permitting much larger step sizes, thus allowing 
solutions to be obtained for problems where the small steps required by the 
explicit technique prevented even the consideration of the case. Finally it 
should be recalled that an extremely simple case was chosen only for purposes 
of illustration and the general technique coded in the modified Program SEAGULL 
will handle many species with multiple reactions. 
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2.5 BOUNDARY CONDITIONS 


By combining the two momentum equations we obtain 



2 P r = 0 


(25) 


where r = u/w . 


The global continuity equation is used in the form 

2 

P + —— (rA . P +T ) - f(s) = 0 . 
z j3 1 r r' A ? (3 


(26) 


Now, using Eqs. (28) and (29) that involve partials of P and j only, a char- 
acteristics compatibility equation can be written 


dr T jr 

dz 2 . 

w ip 


J a i ( ij2 + 

1 T/> 


2 . 
w ) 


dp 

- 1 = RHS 

dz 


(27) 


where the right-hand side (RHS) is: 


RHS = 

A 2 


T 

pijjw \' — 


/ /a i (u 2 +w 2 ) 

t7^ 


- i 


(28) 


Then by combining these equations with the chain rule definitions. 


dP 

dP 

. x 9P 

dz " 

dz 

+ * 3r 

dr 

3 r 

, . 3r 

dz 

3z 

+ * 3? 


(29) 


an equation for 
quantities. 


the axial pressure variation is obtained in terms of known 


P = - X P r + 


w (RHS - T z - X r r ) 


t*/ a i (u 

4J f T /i|i 


2 . 
w ) 


- 1 


(30) 
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where 


, w 2 tAi _ T JAj^ + w 2 ) 

|3 p»pf T7^ 


( 31 ) 


By integrating Eq. (30), the pressure is obtained at the new wall point 
and the velocity is computed from the streamwise momentum equation 


P q 


9q 0P 

9s _ 0S 


(32) 


Velocity components u and w are then obtained from the known wall slope. 
The boundary conditions are completed by integrating the species continuity 
equations along the wall in a manner similar to that used for the interior 
points. 


2.6 SHOCK WAVES 

The general scheme outlined in Appendix A is utilized for computing shock 
wave information. The original ideal gas analysis was supplemented by the 
addition of methodology to handle reacting gas calculations across the shocks. 


Figure 1 illustrates a stream tube passing through an oblique shock 
wave. This wave, which is extremely thin, will cause an almost instantan- 
eous rise in pressure and temperature. For some distance downstream of 
the shock wave (in a reacting gas) a non-equilibrium zone will exist followed 
by a return to chemical equilibrium. The following analysis discusses the 
fluid flow properties in such a way that the non-equilibrium process need not 
be specified in order to arrive at an exact solution for the gas properties. 

It is impossible to determine the location of the new equilibrium shock point 
location without a detailed description of the non- equilibrium reaction process. 
It will be assumed therefore that this zone is thin and that no significant errors 
are introduced by letting the downstream physical location lie on the upstream 
location. 


Consider a control surface as shown in the figure. The conservation 
of mass yields: 

P 2 ^2^2 ~ ^l^l^l = ^ (33) 
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Let ns and ts be unit vectors normal and tangent to the shock surface, the 

unit vectors parallel and perpendicular to the streamline upstream of the 
shock are 


i 0 = sin€ ns + cos£ ts 

i n = -cos€ ns + sin€ ts (34) 



and the unit vectors parallel and perpendicular to the streamline downstream 
of the shock are: 

f 0 ( = sin (6 -6 ) ns + cos (e -6 ) ts 


In' = -cos (6 -6 ) ns + sin (6 -6 ) ts 


(35) 


Conservation of momentum gives 


Pi A i _ 


P 2 A 2 in 


-(P 1 + P 1 q 1 )A 1 lfl+ I5sr in + (p 2 + p 2 q 2 )A 2 ifl- - -~f rr§) 


(36) 


and, after substitution of Eqs. (34) and (35) and setting each component to 
zero, Eq. (36) becomes, 


2 2 
P 2 q 2 cos (6 -8) A 2 - q x 

2 P 1 A 1 
(Pl + Pi q x ) Aj sin« + cos - 


Aj cost = 0 
2 

(p 2 + p 2 q 2 ) A 2 sin (£ -8) 


(37) 


P 2 A 2 COS (£ - S) 
tan (6 - 0) 


(38) 


But from geometry it can be seen that 


tl _ sin(£ -g) 
A^ sine 


(39) 


After substitution of Eqs. (40); (33), (37) and (38) become 


p 2 q 2 sin (€ - 8) - p ^ qj sine = 0 


(40) 
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0 


(41) 


2 2 
P 2 q 2 sin(e-$) cos (6- 5) - p^^ sintcost = 


, 2 . 2 . 2.2 

p 2 + P 2 q 2 Sin < € “ 5) ~ P 1 " P 1 q l sm € = 0 


(42) 


The above set of relations contains €, §, p^, P 2’ q 2 aS un ^ nown quantities, but 


— P ( S 2’ ^2^ J p 2 — P ^ S 2’ q 2^ 


(43) 


So that if one variable, say €, is taken as an independent parameter the 
remaining unknowns (§, q~,s ) may be found by an iterative solution. These 

c 2. 

equations are, of course, formally the same as the ideal gas solution. The 
difference lies only in the variation of pressure etc., with entropy and velocity. 


2.6.1 Iterative Solution of the Oblique Shock Relations 
Rearranging Eq. (44) yields; 


sin(e -$) 


P 1 q 1 sint 

t>Z q 2 


while squaring both sides of Eq. (45) and substituting the above relations 
yields, after simplification; 


q, - q 


2 2 
sin 6 + cos e 


1/2 


= 0 


(44) 


and Eq. (42) becomes 


P 2 + P 1 q l sin f 


- 1 


- P 


= 0 


(45) 
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In functional form Eqs. (44) and (45) are just 


From calculus 


G 1 ( S 2’ q 2^ " 0 
G 2 ^ s 2’ q 2^ = 0 

9G, dG 

dG l = ~dq~ z dq 2 + as^ ds 2 

dG- dG 

dG- — — dq + — ds 0 

2 dq z 2 ds^ 2 


( 46 ) 


(47) 


Now 



(48) 



ds . 


P 


d 

2ds z 




(inp 2 ) 


/ 
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Rather than calculate the partial derivatives numerically by perturbing the 
functions inp^, in,/> approximate values for these derivaties will be found 
by assuming that locally the gas behaves ideally, that is to say 




er 2 <?r 2 


ds~ 5^2 5^2 


= 0 


so that holding q^ constant 


d(ln p 2 ) 

ds„ 


d(fn p 2 ) 

ds 0 


l ~V ; 


in 


and holding s^ constant 


eUnp 2 ) _ 

— ^ it- (to p?) 
dq 2 2 dq 2 2' 


q 2 P 2 

?2 


y 4 c 

writing Eq. (47) in finite difference form: 


( n ) / 

G (n+1) _ G (n) = 

1 1 dq 7 


(n+1) ( 

- q 2 


fs r (^ ) 

"’) * “V- (■ 


(n+1) „ (n) 


- s. 


Ftr ( n ) i \ ( n ) 


Since the root G^ = = 0 is desired, G^ n+ ^\ G 2 ^ n+ ^ are set to zero, 

resulting in 


^where n here is an iteration counter. 


(49) 


(50) 
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(n+1) _ (n) 


= s v “' + G, 


. , dG ■, 

(n) 1 


(n) 




- G 




dG. 


(n) 


dq : 


f dG^ dG^ 


dG ^M dG ^ n ) 


ds 


aq^ 


as. 


aq, 


(51) 


and 




. . 9G 

G, n + -jr- 

i as 


(n j 

1 / s (n+ 1 )_ s ( n) ] 


aG 


(n) 


aq. 


(52) 


The iterative solution using Eqs. (51) and (52) is continued until the desired 
convergence of G^ and G^ is reached. The solution is completed by 


S = € - sin 


p i q i 


p 2 q 2 


sine 


(53) 


The first guess to start the solution is an ideal gas solution to the set of 
equations. If it is indeed an ideal gas under analysis the first guess is exact. 
These relations are 


5 = € - tan 


-1 


tane 


1 


2.2 

\M l sine 


+ 


V 1 


V 1 


q? = 


cose 

q l cos (€-8) 


(54) 


s 2 = + 


R 1 

in 

2X M 2 sin 2 - (Y - 1) 

+ X x in 

- 

tan(e-S) 


V 1 ' 

y i + 1 

tan€ 





_ _ 



20 



2.7 CONTACT SURFACES 


Contact surfaces or slip lines are computed utilizing the condition 
that the pressures and flow angles must match across the surface. The 
equations for pressure and flow angle developed for the wall boundary con- 
ditions are also applied at the contact surface except that the streamline 
slope must now be obtained iteratively utilizing the matching conditions. 
Equations (25) and (26) are transformed to the computational coordinates 
and integrated via the MacCormack operator as described in Appendix A. 

The predictor (first level) is applied to the high Mach number side of the slip 
line and the corrector (second level) is applied to the low Mach number side. 

With the pressure and slope now known at the contact surface the 
velocity is again obtained from Eq. (35). The velocity components are then 
evaluated from the known slip line slope. 

2.8 THERMODYNAMIC PROPERTIES 

The modified Program SEAGULL handles general gases composed of 
various chemical species rather than a global gas. This requires thermo- 
dynamic data inputs for all participating species as well as for all vibra- 
tional levels of the various excited molecules involved in vibrational energy 
exchange. For these molecules it is assumed that the molecular rotation 
remains in equilibrium with the translational temperature but each molecule 
is allowed to vibrate independently, corresponding to the energy in its re- 
spective vibrational level. 

Thermodynamic data are supplied to the code in the form of tabular 
inputs of specific heat, entropy and enthalpy for each species as a function 
of temperature and are taken directly from the JANAF thermochemical 
tables (Ref. 4). 
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2.9 VIBRA TIONA L NONEQUILIBRIUM 


/ 


To make maximum use of the current technology of finite rate formula- 
tions, the various vibrational levels of the respective excited molecules are 
treated as individual chemical species. In order to treat them as individual 
species, it is first assumed that the excited molecule can be resolved into its 
various vibrational levels each of which is then treated as a separate species. 
This assumption requires that the thermodynamic properties be separately 
specified for each vibrational level. It is further assumed that molecular 
rotation remains in equilibrium with the translational temperature but the 
molecules can vibrate independently of the translational temperature. Ac- 
cordingly, each molecule is permitted to vibrate independently corresponding 
to the energy in its respective vibrational level. 
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Appendix A 


SHOCK FITTING METHOD FOR COMPLICATED 
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Shock Fitting Method for Complicated 
Two-Dimensional Supersonic Flows 

Manuel D. Salas* 

NASA Langley Research Center, Hampton, Va. 


The floating shock fitting technique is examined. Second-order difference formulas are developed for the 
computation of discontinuities. A procedure is developed to compute mesh points that are crossed by discon- 
tinuities. The technique is applied to the calculation of internal two-dimensional flows with arbitrary number of 
shock waves and contact surfaces. A new procedure, based on the coalescence of characteristics, is developed to 
detect the formation of shock waves. Results are presented to validate and demonstrate the versatility of the 
technique. 


I. Introduction 

I NVISCID supersonic flows are governed by a first- 
order quasi-linear hyperbolic system of equations. The 
numerical computation of these flows is complicated, because 
the regular solution may break down due to the nonlinearity 
of the governing equations. The breakdown is characterized 
by the appearance of surfaces, such as shock waves and vortex 
sheets, across which the dependent variables or their 
derivatives are discontinuous. During the last decade, 2 
numerical techniques have emerged for the analysis of these 
flows. One, known as shock capturing, tries to remove the ex- 
plicit computation of the discontinuties by generalizing the 
concept of a solution of the Euler equations to include weak 
solutions (i.e. , discontinuities). Because the shock capturing 
scheme requires no special treatment to deal with discon- 
tinuties, it has become an extremely popular way of com- 
puting. However, despite its present popularity, the results 
obtained with this technique 1 force us to agree with Moretti’s 
conclusion that shock capturing is a poor interpretation of a 
physical phenomenon, and an extremely uneconomical way of 
computing . 3 

The other technique, known as shock fitting, makes special 
provisions for explicitly computing the discontinuties. In 
essence, it locates the discontinuties and treats them as boun- 
daries between regions where a regular solution is valid. The 
effectiveness and soundness of this approach has been proven 
for a large number of problems involving 2, 3, and 4 in- 
dependent variables. However, for problems in 3 or more in- 
dependent variables, the partitioning of a flowfield into 
regions where a regular solution is valid, can create some dif- 
ficult topological problems . 3 Recently, Moretti 4 has 
developed a technique (reminiscent of that proposed by Richt- 
myer and Morton , 5 that treats the discontinuities explicitly 
but does not require them as boundaries of the flow. This 
therefore eliminates the problems associated with partitioning 
the flowfield. It is this technique, known as floating shock fit- 
ting, that forms the subject matter of this paper. 

To avoid unnecessary complications, the technique will be 
discussed in the context of a two-dimensional problem. 
However, to demonstrate the capabilities of the technique, 
results from complex flowfields are presented. 


II. Problem Definition 

We will discuss the floating shock fitting technique as it ap- 
plies to the flow of an inviscid, perfect gas through a two- 
dimensional duct of arbitrary geometry. The flow is assumed 
to be supersonic throughout, and consequently, irregular 
shock reflections are not considered. Shock waves generated 
by either discontinuties in the slopes of the walls, or by 
coalescence of pressure waves, will be automatically detected 
and fined into the calculation. The number of shock waves, 
contact surfaces, and the number of interactions of these 
discontinuties, is considered arbitrary (limited only by com- 
puter storage). 

III. Governing Equations 

The flow variables arc nondimensionalized by scaling the 
pressure, density, and temperature, with respect to their 
initial values (p„ p„ T,) , scaling the velocities with respect to 
( p , /p,) and scaling all lengths with respect to some ar- 

bitrary reference length. The Euler equations written in a Car- 
tesian frame (r,z), with thez-axis running along the length of 
the duct, are 

P- + ( uwP r + ywu r — y wvv r ) //3 = 0 
tv. + (uww r — yTu r - TuP r ) /(3 = 0 
u. + ( uu , + TP r ) /w — 0 

S.+uS r /w — 0 (1) 

Where P is the natural logarithm of pressure, u and vv are the 
velocity components in ther and z directions, respectively; S is 
the entropy, 7 is the ratio of specific heats; the temperature is 
given by 

7"=exp((S+ (y — J) P) /y) (2) 

and 0 = w 2 — yT. 

The following equations are introduced here for later 
reference. By combining the 2 momentum equations, we get 
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7 r + 


W 2 TT r 

0 


+ (/ + 


VA=0 


(3a) 


where 7 = « /w. With Eq. (3a) and the continuity equation in 
the form 


P.+ w 2 (rP r + yT r )/0 = O (3b) 
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we can obtain the compatibility equation 


t' ± ( [yT((5+ w 2 t 2 ) ] ,/: /yw 2 )P' =0 (4) 

where the primes denote differentiation along the charac- 

teristic directions 

X* = (w 2 r±[yTU3+w 2 T 2 )} *)//3 (5) 

The Euler Eqs. (1) are recast in terms of computational 
coordinates (X.Z) defined by 

X={r-b{z))/(c(z)-b(z)) Z=z (6) 

where r=b(z) and r=c(z) define the lower and upper walls 
of the duct. In matrix form the equations of motion become 

fz=4fx (7) 


h — az— H 



; 

A^-S. 
N +2 

i 1 
AX 
♦ , 


nTi 

> k) 



,N - 1 


. j 

,N - 2 

NNWVXWWWVvX 

K K+l 


Fig. 1 Computational plane. 



Fig. 2 Computation in the neigh- 
borhood of a shock wave. 


where the vector/ is 


level 2 


r p 


fhl,K+l “ 1/2 if N.K +/n,K + 1 +4n.K + 1 


f= 


u 

w 

S 


and the matrix A is 

I - a n a, 2 a u 0 


llK+l )^Z/AX) (10) 

Because the discontinuties are allowed to float over the 
fixed grid, special provisions must be taken when evaluating 
the X derivatives of /, at mesh points in the neighborhood of a 
discontinuity. Consider again the computation of point 
(N.K+l ) , but now in the presence of a shock as shown in Fig. 
2. There is no difficulty in computing the first level of the 
MacCormack scheme, since the forward difference 


«2! 

an 

0 


a 22 

a 32 

0 


0 

a,, 

0 


0 

0 


a n ~ — (uwX r + 0X Z ) /0 
a 12 = —ywX r /@ 


(9) 


f x — Wn + i.k ~/n.k) / ' &X 

is valid. However, for the second level, the backward dif- 
ference 


fx — if n.k+i ~/n-i,k+i) f 

is not allowed because /is not differentiable across the shock. 
On the other hand, the very simple approximation, 


a u =yuX r /f3 


fx - (/n.k+i f s,K + 1 ) / (X N —X sK+ ,) 


a n = - TX r /w 
a 22 = — (uX r + wX . ) /w 
a 3 , =uTX r /(3 
a 32 =yTX r /(S 

IV. Interior Point Computation 

Unlike shock fitting, floating shock fitting does not sub- 
divide the computational plane into regions bounded by 
discontinuities; instead, a single region is used, as illustrated 
in Fig. 1 . In this single region, mesh points are equally spaced, 
and all discontinuties are allowed to move freely or float over 
the fixed grid. Their motion is subject only to the governing 
equations. 

Except for the mesh points at the walls and in the neigh- 
borhood of a discontinuity, all mesh points are computed 
using the MacCormack scheme 6 to integrate Eq. (7). For 
example, consider the section of a duct, in the computational 
plane, sketched in Fig. 1. Using the MacCormack scheme, 
and assuming the solution is known at the ATth steps, /*,*+/ is 
obtained from 

level 1 

/n,k+ i = f n,k + 4n.k if i.k ~/n,k ) h2.f 


where s denotes conditions at the shock, can result in large 
truncation errors, and the computation becomes unstable. 
Moretti recommends using an approximation that involves 
the values at the shock and the 3 neighboring mesh points 
[Eq. (61) in Ref. 2] . We have not used this approximation, 
because it makes the numerical domain of dependence point 
( N,K+1 ) much larger than it physically is. Instead, we use 
the following approximation 

fx = ~ f (0.5e— l)fs+i'K+i +0-5f s ,K + i 

+ 0.5 (1 — e)/ N K + / 1 / AX (II) 

where 

e=iX NK + l -X,_ K+l )/AX 

The truncation error for Eq. (1 1), (based on f x , not f) is given 
by 

(e + 2)(e-l)AX/4f xx 

With Eq. (11), the computation of point (N,K+1) proceeds 
as the computation of any other interior point. Similar con- 
siderations apply to the computation of point ( N — I.K + 1). 

If point ( N+I , K+l) is not available when applying Eq. 
(11), as is the case indicated in Fig. 3 where 2 shocks are 


25 


Fig. 3 Computation with 
closely spaced shocks. 



Here, the X derivatives are evaluated as three-point one-sided 
differences and 




b zt for X" 
Cf. for X + 


Then, following the MacCormack algorithm, the value of the 
pressure at the wall, P w , is given by 


Fig. 4 Compulation when a shock 
crosses a mesh line. 



shown close to each other, then if the distance between the 
shocks is greater than AX, / x is approximated by 

f X — if: s2. K + I ~fsl.K+l ) 

+ i — X SI K + t ) (12) 

If, however, the distance between the shocks is smaller than 
AX, the values at (N.K + l) are interpolated from the shock 
values. 

When a discontinuity crosses a mesh line, the computation 
of this mesh point must be modified. In Fig. 4, an up-running 
shock is shown crossing the Mh mesh line between steps K 
and K+l. The usual interior point computation would 
require taking the one-sided difference between points (N+K, 
K) and ( N,K) for level one. This procedure is no longer 
proper, because of the presence of the shock. We, therefore, 
proceed by evaluating points A and B as follows 

f A ~ P N— 1,K T ^ C/*/V— l.K = 1 f N— /.X ) 1 

f B —fs.K + 5 ( fs,K+ I f s,K ) J (13) 

where 

5=(X N -X^)/(X S ' K + l -X, iK ) (14) 

then, the value of /for level one at (N, K+l) is 

f\.K + i = f b + dfi C/o f a ) ( I — 5) AZ / AX (15) 

Care must be taken in evaluating 5 when dA\/dZ = 0., since 
Eq. (14) becomes indeterminate. We therefore require that, 


P »-.A C+/—P it'. K+P Y.&Z 

P»,K+, = 'A{P^ + P^ + l +Py_AZ) 

With the pressure and entropy known at the wall, the tem- 
perature is obtained from Eq. (2), and the magnitude of the 
velocity from 


u‘ + W 1 7 

4- T= constant (17) 

2 y-J 

The two-velocity component can then be determined from the 
known wall slope and the vanishing of the velocity component 
normal to the wall. 


VI. Shock Detection 

Perhaps the most distinguishing feature of the problem we 
are considering, and certainly its most troublesome one, is the 
possible formation of shock waves. In this section, we will 
consider the procedure required to detect an incipient shock. 
There are two important features for any such procedure. 
First, the procedure must be reliable. Incipient shocks must be 
located at the proper place, but most important, the 
procedure must be able to interpret the physics of the problem 
accurately, and should not be so sensitive as to create 
nonexistent shocks. Secondly, because the procedure must be 
executed at every step, it must be very fast and require a 
minimum amount of data to make a reliable decision. 

The procedure we have constructed, based on the concept 
of coalescing characteristics, requires only information on 2 
adjacent mesh points, and we find it to be both reliable and 
fast. However, its empirical nature should be kept in mind. It 
is applied at the end of a computational step as follows: 

Assume that the flowfield is known at the Xtth step (see Fig. 
5). At each mesh point /V and /V— / , the up-running (Xf,X.* 
v _ / ) and downrunning (X^, X ) characteristics are 
evaluated using Eq. (5). If we now assume that the AZ in- 
tervals between steps (evaluated from the CFL condition)' 1 are 
essentially the same, then a quantity rj±, which is the 
reciprocal of the number of steps needed for 2 characteristics 
to coalesce can be evaluated from 


i dA\ | 


AZ > 0.2 AX 


' dZ 

otherwise, 6 is set to one to give the desired result 
J N,K + ! = fs.K+l 

A similar computation is performed when a down-running 
discontinuity crosses a mesh line. 


Pz=~dP x 

where 


V. Wall Computation 

>n, E 

yw 2 


The compatibility equation, Eq. (4), is applied at the wall to 
evaluate P z 


[yT($+ W 2 r 2 )} 

X=A'.-+X ± A r , 


2 , , ./, (t z + \t x ) (16) 


t) ± — ( X£_/ — X v ) Ac/ A r 


A considerable number of numerical experiments indicate 
that a shock should be fitted into the calculation when 17 * 
becomes greater than 0.22 (approximately 4 steps before the 
characteristics coalesce). 

We have found that for a wide range of problems, if the 
detection of incipient shocks is delayed until jj* < 0.22), (hen 
the procedure becomes very sensitive and shocks are predicted 
unnecessarily. 


Fig. 5 Details of the shock 
detection. 
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Fig. 6 Evaluation 
of (he low-pressure 
side of a shock. 
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The detected shock is assumed to have the average slope of 
the coalescing characteristics and to be located half way bet- 
ween r N and r N _,. With the initial shock location known, the 
initial flowfield variables on the low pressure side are linearly 
interpolated between the 2 mesh points bracketing the shock, 
and on the high pressure side they are linearly extrapolated 
using the 2 mesh points next to the shock on the high pressure 
side. 

Shocks that are generated at the walls of the duct, because 
of slope discontinuities at the walls, are introduced from user 
supplied information on the wall geometries, and do not make 
use of the procedure described in this section. 


VII. Shock Computation 

The low-pressure side of a shock is computed using the 
same scheme described for the interior point computation, but 
using the X derivatives given in this section. In Fig. 6, 
assuming that mesh point N lies on the low-pressure side, the 
X derivatives at the shock for the predictor level are evaluated 
from 


fx s -[ 


2 

7+7 


(A* Sn - /.A ) +0.5/n - 2, K 


(0.5 + e) 
U+2f) 



/AX 


( 18 ) 


where 


t={X^ K -X N ' K )/AX (19) 

If only two mesh points are available, then the derivatives are 
evaluated from 


fx s ~ ["/ + f /tv -/.*■+ ( A. k ~ 2f n . a j / &X (20) 

with t given by Eq. (19). The truncation error for Eq. (18) is 
of order AX 2 and for Eq. (20) it is 


(e 


n 



It should be noted that Eqs. (1 8-20) are not the same formulas 
suggested by Moretti, 2 ' 4 particularly Eq. (18) has a lower 
truncation error than Moretti’s counterpart. 

If only 1 mesh point is available (Fig. 6b), then f Xs is ap- 
proximated by 


fx s ~ C/l/.A ~fs2.K ) / (X s /A — X s2 .k ) (21) 

provided that (X sUi ~X, 2K ) is greater than AX. If (X, LK - 
X s2K ) is smaller than AX, then f Xs is set to zero. For the 
corrector stage, the same derivatives are used but with / 
replaced by f. 

If the low-pressure side lies on the other side of the shock, 
the derivatives that are used are the negatives of those given 
by Eq. (18) and Eq. (20), and c is given by the negative of Eq. 
(19). 


The computation of the high-pressure side of the shock has 
been formulated using the compatibility equation, Eq. (4), on 
the characteristic reaching the shock on the high-pressure side 
along with Rankine-Hugoniot conditions. The details of this 
classical treatment of the shock can be found in Ref. 4. 

The computation of the high-pressure side was the only 
computation requiring involved programing logic. The main 
difficulty is the determination of the origin of the charac- 
teristic reaching the shock in the immediate vicinity of other 
discontinuities. It should be, however, possible to avoid these 
difficulties by computing the high-pressure side with the 
method suggested by Kentzer. 7 

VIII. Contact Surface Computation! 

To compute a contact surface, we exploit the fact that both 
pressure and streamline slope are continuous across a contact 
surface. Eqs. (3a, b), which involve derivatives of P and t 
only, are thus recast in terms of the computational coor- 
dinates ( X,Z ), and are integrated using the MacCormack 
scheme where the X derivatives are evaluated using Eq. (20). 
In the predictor level, Eq. (20) is applied on the high Mach 
number side of the contact surface. In the corrector level, it is 
applied on the low Mach number side. This choice is ar- 
bitrary, and the same results are obtained if the order is re- 
versed. 

With the pressure and entropy known at the contact sur- 
face, the magnitude of the velocity follows from Eq. (17), and 
the 2 velocity components can be determined from the 
evaluated streamline slope. 

IX. Interactions 

The interactions of the discontinuities with the walls of the 
duct and with each other are computed by locally exact 
solutions. The details of these calculations can be found in 
most text books in fluid mechanics, for example. Ref. 8, and 
will not be discussed here. 

X. Elimination of Weak Discontinuity 

Discontinuities with a Mach number difference across the 
discontinity of less than 1% are automatically eliminated 
from the computation. This feature is not a requirement of 
the technique. It has, however, been incorporated, because 
for most practical applications, tracing these weak discon- 
tinuties is not necessary. It speeds up the computation, and it 
avoids computing degenerate discontinuities (i.e. the contact 
surface formed by the interaction of 2 shock waves of equal 
strength). 


XI. Results 

The results of several calculations performed with the 
floating shock fitting technique are presented to demonstrate 
the accuracy and versatility of the technique. All the com- 
putations were done on a CDC 6600 computer. The core 
required to run the code is under 55000 8 locations. 

The flowfield produced by a uniform flow impinging on a 
wedge is often used to validate numerical calculations. The 
problem is a trivial one, consisting of uniform properties in 
front and behind the shock wave. The exact solution for a 
freestream Mach number of 3.50 which is initially compressed 
by a 15° wedge and then expanded, such as to exactly cancel 
the shock wave generated by the wedge, is shown in Fig. 7. In 
the figure, the lower wall pressure distribution shown 
corresponds to the case where the expansion occurs a very 
short distance before the shock reflects from the wall. The 
exact solution is obtained with the floating shock fitting 


fWe have recently experienced some difficulties with this procedure 
when the contact surface is next to a strong expansion fan. The 
problems have been overcome by reformulating the contact surface 
compulation along the guide lines discussed in the shock compulal ion. 
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Fig. 7 Flowfield over a 
15° wedge wilh shock 
cancellation. 



Fig. 10 Upper wall 
pressure for M = 2 
case. 



Fig. 8 Parametric 
study of a nozzle, 
showing shock 
waves and isobars. 



Fig. 9 Shock path for 
M= 1 case. 



technique using only 2 mesh points in less than 1 sec of com- 
putational time (the same results are obtained increasing the 
number of mesh points). 

In Fig. 8, the results of a parametric study of the flowfield 
produced by a nozzle are presented. The initial Mach number 
was used as the parameter, allowing it to range from near 1 to 
4. Both the isobars and shock-wave patterns are shown in the 
figure. For the Mach 4 case, most of the isobars downstream 
of the shock wave were not included in the figure to avoid 
overcrowding of lines. The total computational time for these 
4 cases, using 30 mesh points, was approximately 25 sec. 
Figures 9 and 10 give some further details in the calculation of 
the Mach 2 case. In these 2 figures the results obtained with 5, 
10, 20, 30, and 60 points are shown. Figure 9 shows the 
calculated shock-wave path, while Fig. 10 shows the upper 
wall pressure distribution. The results show that a good 
qualitative result for this more complicated problem can be 
obtained with as few as 5 mesh points, and that good quan- 
titative results arc obtained with only 1 0 mesh points. 

It is very often argued that shock fitting is limited in 
its capability to handle flowfields with more than a single 
shock. 911 In Refs. 3 and 12 this has been shown not to be the 
case for very complex three-dimensional steady, and one- 
dimensional time-dependent flows. We hope to fill the 
existing gap for two-dimensional steady problems with Fig. 1 1 
and also to demonstrate that floating shock fitting is a 
numerically stable technique even in multishocked flowfields. 
In Fig. 11, the geometry of NASA Langley’s scramjet com- 
bustor model is simulated. The results shown 
are for a flow with an incoming Mach number equal to 3. 
inclined at 8° with the horizontal. The resulting flowfield con- 



Fig. II Flowfield for a simulated .scramjet, showing shock waves, 
vortex sheets and isobars. 



Fig. 12 Plume flowfield showing shock waves and vortex sheets. 

sists of a complex pattern of interacting expansions and shock 
waves. In the figure, the shock-wave pattern is shown by the 
heavy lines. The contact surfaces are shown by the dashed 
lines, and isobars are shown by the light lines. Certain discon- 
tinuities are seen to terminate abruptly in the flowfield, while 
some shock interactions do not show the expected vortex 
sheet; this is a result of the automatic elimination of weak 
discontinuties previously discussed. The calculation was per- 
formed in 6 pieces, using a total of 120 mesh points and 
requiring approximately 2.2 min of computational time. The 
calculation shows approximately a dozen discontinuties, with 
some 40 interactions. 

The technique has been extended to internal axisymmetric 
flows (with and without swirl), internal line source type flows, 
and to external two-dimensional flows. An example of the ap- 
plication to external flows is shown in Fig. 12, where a flow 
with an initial Mach number of 2.0 exhausts into a Mach 3.0 
freestream. The lower wall ends at z = 4.0, while the upper 
wall ends at z = 4.2. The freestream conditions for the ex- 
ternal flows are shown in the figure. 

XII. Conclusions 

Floating shock fitting for two-dimensional steady flows is 
an accurate and stable technique capable of solving very'Com- 
plex problems with acceptable running times. A fully three- 
dimensional code is planned for the near future. 
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INPUT GUIDE AND SAMPLE CASES 


DESCRIPTION OF INPUT DATA 

ALL INPUT DATA IS READ FROM TAPE 5. 

*THE FIRST INPUT CARD CONSISTS OF 11 INTEGERS. 

THE FORMAT IS 1615 

THESE INTEGERS ARE ID, KOUT, ICN, IBN, ISO, NA, ISY, JJ, ITY, ITP, IVI 
ITY IS USED TO INDICATE MERGING OF INTERNAL DUCTS, AND IS NOT 
RELATED TO JETS OR PLUMES. 

ID=0 INPUT CARD NO. 3 IS NOT READ 
ID= 1 INPUT CARD NO. 3 IS READ 

KOUT = THE NO. OF STEPS TAKEN BEFORE A COMPLETE OUTPUT OF 
THE FLOW FIELD IS MADE 

ICN = THE TOTAL NO. OF SHARP CORNERS ON THE UPPER WALL 
ICN CANNOT EXCEED 10 

IBN = THE TOTAL NO. OF SHARP CORNERS ON THE LOWER WALL 
IBN CANNOT EXCEED 10 
ISO=0 NO ISOBAR OUTPUT WILL BE PROVIDED 

ISO=l ISOBARS WILL BE COMPUTED AND PRINTED AT EVERY STEP 
NA = INITIAL NUMBER OF MESH POINTS MINUS TWO 
NA CANNOT EXCEED 98 

IN MOST CASES 20 POINTS ARE SUFFICIENT, AND ON THE AVERAGE 
THERE SHOULD BE ABOUT 5 POINTS BETWEEN DISCONTINUITIES. 

ISY=0 NO SYMMETRY IS ASSUMED IN THE COMPUTATION OF THE LOWER 
WALL 

ISY= 1 SYMMETRY IS IMPOSED IN THE COMPUTATION OF THE LOWER WALL 

JJ=0 FOR TWO DIMENSIONAL FLOW 

JJ=1 FOR AXISYMMETRIC FLOW 

JJ= 2 FOR LINE SOURCE FLOW 

ITY=0 A SINGLE DUCT CALCULATION 

ITY= 1 SINGLE DUCT CALCULATION WITH DATA WRITTEN ON TAPE 7 
FOR FUTURE MERGING RUN. CAN ALSO BE USED TO RESTART 
A CALCULATION 

ITY=2 DATA IS READ TO MERGE TWO OR MORE FLOWS, AND DATA IS 
WRITTEN ON TAPE 7 FOR FUTURE CONTINUATION 
ITY=3 DATA IS READ TO MERGE TWO OR MORE FLOWS 
ITP=0 NO CALCOMP PLOT IS GENERATED 
ITP= 1 A CALCOMP PLOT IS GENERATED 
IVLO NO VISCOUS FORCES ARE COMPUTED 

IVI=1 SPALDING -CHI USED TO COMPUTE VISCOUS FORCES, THE 
REYNOLDS NO. WILL BE BASED ON THE WALL GEOMETRY 
LENGTH 
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*INPUT CARD NO. 2 CONSISTS OF 5 INTEGERS 
THE FORMAT IS 1615 
NS = NO. OF CHEMICAL SPECIES 

NT = NO. OF POINTS IN THERMODYNAMIC TABLES 
NM = NO. OF CATALYTIC SPECIES 
NR = NO. OF REACTIONS IN MECHANISM 
IFROZE=0, FINITE RATE 
IFROZE= 1 , FROZEN FLOW 

*INPUT CARD NO. 3 IS READ ONLY IF ID= 1 
ITS PURPOSE IS SOLELY FOR IDENTIFICATION, A FREE FIELD FORMAT 
IS USED BETWEEN COLUMNS 6 THRU 55 

*INPUT CARD NO. 4 CONSISTS OF 5 FLOATING POINT NUMBERS 
THE FORMAT IS 8E10.4 

THESE NUMBERS ARE ACH, GM, DIST, TOLER, SCALE 

ACH = INITIAL MACH NUMBER 
GM = RATIO OF SPECIFIC HEATS 

DIST = DISTANCE TO THE END OF THE CALCULATION (CM) 

TOLER = A TOLERANCE PARAMETER TO INCREASE OR DECREASE THE 
NO. OF MESH POINTS. IT WORKS AS FOLLOWS: IF SET TO X 
THEN THE NO. OF MESH POINTS WILL BE INCREASED WHEN 
THE DISTANCE BETWEEN MESH POINTS IS X TIMES THE INITIAL 
DISTANCE BETWEEN MESH POINTS, AND MESH POINTS ARE DE- 
CREASED WHEN THE DISTANCE BETWEEN MESH POINTS IS l./X 
TIMES THE INITIAL DISTANCE BETWEEN MESH POINTS. 

SETTING TOLER TO A LARGE NUMBER (SAY 1000.) HAS THE 
EFFECT OF NQT CHANGING THE NUMBER OF MESH POINTS. 
SCALE = IS A SCALE FACTOR FOR CALCOMP PLOTS NEEDED ONLY IF 
ITP= 1 BOTH THE R AND Z COORDINATES ARE MULTIPLIED 
BY SCALE BEFORE THEY ARE PLOTTED. 

-INPUT CARD NO. 5 CONSISTS OF 7 FLOATING POINT NUMBERS 
THE FORMAT IS 8E10.4 

THE NUMBERS ARE ZEXIT, PINF, TINF, GMINF, ACHINF, ANGINF, RCINF 
THIS INFORMATION IS NEEDED IN CASE THE LOWER WALL OF THE DUCT 
OPENS TO THE FREE STREAM AT Z = ZEXIT. IF THIS IS NOT THE CASE 
THEN ZEXIT SHOULD BE SET TO A NUMBER GREATER THAN THE VALUE 
OF DIST, AND THE REMAINING SIX NUMBERS NEED NOT BE DEFINED. 

ZEXIT = LOCATION OF EXIT PLANE (CM) 

PINF = FREE STREAM STATIC PRESSURE (ATM) 

TINF = FREE STREAM STATIC TEMPERATURE ( K) 

GMINF = NOT REQUIRED PRESENTLY 
ACHINF = FREE STREAM MACH NUMBER 

ANGINF = FREE STREAM INCLINATION IN DEGREES, MEASURED POSITIVE 
COUNTERCLOCKWISE FROM THE Z AXIS 
RCINF = NOT REQUIRED PRESENTLY 


♦INPUT CARD NO. 6 (REQUIRED ONLY IF ZEXIT. LE. DIST) 

THE FORMAT IS 8E10.4 

THIS CARD SUPPLIES THE SPECIES DISTRIBUTION IN THE EXTERNAL 
STREAM IN THE SAME ORDER AS IN THE MAIN STREAM* 

ALFA(1, I), FIRST SPECIES IN THERMO DATA TABLES 
ALFA(2, 1), SECOND SPECIES IN THERMO DATA TABLES 
ALFA(N, 1), N-TH SPECIES IN THERMO DATA TABLES 

♦INPUT CARD NO. 7 IS THE SAME AS CARD NO. 5, BUT FOR THE CASE 
WHEN THE UPPER WALL OF THE DUCT OPENS TO FREE STREAM. 

♦INPUT CARD NO. 8 IS THE SAME AS CARD NO. 6 BUT FOR THE CASE 
WHEN THE UPPER WALL OF THE DUCT OPENS TO FREE STREAM. 

♦INPUT CARD NO. 9 CONSISTS OF ICN FLOATING POINT NUMBERS, IF 
ICN=0 THE INPUT CARD IS NOT NEEDED. THESE NUMBERS ARE 
STORED IN ARRAY ZIC(N). 

THE FORMAT IS 8E10.4 

ZIC(N) = LOCATION OF THE NTH CORNER ON THE UPPER WALL (CM) 

♦INPUT CARD NO. 10 CONSISTS OF ICN FLOATING POINT NUMBERS, IF 
ICN=0 THE INPUT CARD IS NOT NEEDED' THESE NUMBERS ARE STORED 
IN ARRAY AIC(N) . 

THE FORMAT IS 8E10.4 

AIC(N) = THE CHANGE IN ANGLE IN DEGREES AT THE NTH UPPER WALL 
CORNER. AIC IS DEFINED POSITIVE IF THE CORNER IS A COM- 
PRESSION AND NEGATIVE IF AN EXPANSION. NOTE THAT IT IS 
A CHANGE IN ANGLE AND NOT THE NEW WALL ANGLE. 

♦INPUT CARD NO. 11 IS THE SAME AS CARD NO. 6, BUT FOR THE LOWER 
WALL. THE DATA IS STORED IN ARRAY ZIB(N), THE INPUT CARD IS NOT 
NEEDED IF IBN=0. 

♦INPUT CARD NO. 12 IS THE SAME AS CARD NO. 7, BUT FOR THE LOWER 
WALL. THE DATA IS STORED IN ARRAY AIB(N), THE INPUT CARD IS NOT 
NEEDED IF IBN=0. A IB IS DEFINED IN THE SAME WAY AS AIC. 

♦INPUT CARD NO. 13 CONSISTS OF 6 FLOATING POINT NUMBERS 
THIS CARD IS NOT NEEDED IF IVI=0 
THE FORMAT IS 8E10.4 

ATWTI(l) = NONDIMENSIONAL LOWER WALL TEMPERATURE 
ATWTI(2) = NONDIMENSIONAL UPPER WALL TEMPERATURE 
XL( 1) = LOWER WALL VIRTUAL ORIGIN FOR SPALDING-CHI (GREATER 
THAN 0.) 

XL(2) = UPPER WALL VIRTUAL ORIGIN FOR SPALDING-CHI (GREATER 
THAN 0.) 

REI = REYNOLDS NO. PER UNIT LENGTH 
PR AND = PRANDTL NO. 


*IF PARTICULAR SPECIES DOES NOT EXIST ENTER ZERO IN THAT LOCATION. 



♦INPUT CARD NO. 14 INITIAL. PRESSURE AND TEMPERATURE IN INTERNAL 
DUCT. 

THE FORMAT IS 8E10.4 

PIN = INITIAL PRESSURE (ATM) 

TIN = INITIAL TEMPERATURE (°K) 

♦INPUT CARD(S) NO. 15 SPECIES DISTRIBUTION MOLE FRACTION IN 
INTERNAL DUCT. 

ALFA{ 1, 1), FIRST SPECIES IN THERMO DATA TABLES 
ALFA(2, 1), SECOND SPECIES IN THERMO DATA TABLES 
ALFA(N, 1), NTH SPECIES IN THERMO DATA TABLES 

♦INPUT CARD(S) NO. 16 SPECIES ID AND THERMODYNAMIC DATA 
THE FORMAT IS A6, 2E10.3 

AID(I), ITH SPECIES IDENTIFICATION 

WTM(I), ITH SPECIES MOLECULAR WEIGHT (g/g-mole) 

HF(I), ITH SPECIES HEAT OF FORMATION (K-cal/mole) 

TTB(J), JTH TEMPERATURE POINT IN TABLES (K) 

CPTB( 1, J) SPECIFIC HEAT FOR ITH SPECIES AT JTH TEMPERATURE 
POINT (cal/ mole -K ) 

STB(K, J) ENTROPY FOR ITH SPECIES AT JTH TEMPERATURE POINT 
(cal/mole-K) 

HT B (I , J ) ENTHALPY FOR ITH SPECIES AT JTH TEMPERATURE POINT 
(K -cal/mole) 

♦INPUT CARD(S) NO. 17 CATALYTIC SPECIES DATA (OPTIONAL) INPUT ONLY 
IF NM=0 

THE FORMAT IS A6/(16E5.2) 

IF NM/O, INPUT THE FOLLOWING 
AID(J), NAME OF JTH CATALYTIC SPECIES 

WF ( J, K), WEIGHTING FACTOR OF ITH SPECIES FOR JTH CATALYTIC SPECIES 
AID(J), NAME OF JTH CATALYTIC SPECIES 

WF(J, I), WEIGHTING FACTOR OF ITH SPECIES FOR JTH CATALYTIC SPECIES 

♦INPUT CARD NO. 18 NAME OF REACTION SET 
THE FORMAT IS 13A6, A2 

♦INPUT CARD NO. 19 REACTION RATE DATA (OPTIONAL) INPUT ONLY IF 
IFROZE IS ZERO 

THE FORMAT IS A6, IX, A6, IX, A6, IX, A6, IX, A6, IX, A6, 7X, 12, II, E8.2, 
F5.2, F10.1, 1F6. 2 
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AID(J), SPECIES (NAME) APPEARING IN ITH REACTION 
IRR(I), REACTION TYPE 

IRT(I), RATE CONSTANT TYPE * 

RCT(I, N), RATE CONSTANT COEFFICIENTS N= 1,4 


FOR MOST RUNS, WITH THE EXCEPTION OF CASES WHERE ITY=2 OR ITY=3, 
THIS IS ALL THE INFORMATION REQUIRED AS INPUT DATA. CASES WHERE 
ITY=2 OR ITY= 3 WILL BE DISCUSSED LATER. LET US CONSIDER NOW THE 
DEFINITION OF THE WALL GEOMETRIES. 


THE LOWER WALL GEOMETRY SHOULD BE DEFINED IN SUBROUTINE 
GEOMB. THE USER DEFINES THE GEOMETRY BY GIVING R OF THE 
WALL AS A FUNCTION OF THE AXIAL COORDINATE Z. THE CODE 
ASSUMES THAT THE INITIAL STATION IS AT Z=0. IF IT IS NECESSARY 
TO START THE CALCULATION AT SOME OTHER VALUE OF Z, THIS CAN 
BE DONE BY REDEFINING Z IN THE MAIN PROGRAM. IN ADDITION THE 
FIRST DERIVATIVE, RZ, AND THE SECOND DERIVATIVE, RZZ, ARE ALSO 
REQUIRED. 

THE UPPER WALL GEOMETRY SHOULD BE DEFINED IN SUBROUTINE 
GEOMC. THE GEOMETRY IS DEFINED IN THE SAME MANNER AS DE- 
SCRIBED ABOVE. A SAMPLE UPPER GEOMETRY FOR A DUCT WITH 
AN INITIAL HEIGHT OF ONE UNIT, A 10. DEGREE COMPRESSION AT 
Z= .5 AND A 6 DEGREE EXPANSION AT Z=l. IS: 

R= I . $ RZ=0. $ RZZ=0. $ IF (Z.LT.0.5) RETURN 

R=l. -. 17633*(Z-.5) $ RZ= -.17633 $. RZZ=0. 

IF (Z.LT.l.) RETURN 

R=.911835+.1051*(Z-1.) $ RZ=.1051 $ RZZ = 0. $ RETURN 

CARE SHOULD BE EXERCISED WHEN DEFINING WALL GEOMETRIES, AS 
THE MOST FREQUENT CAUSE OF PROBLEMS THAT USERS HAVE EX- 
PERIENCED WITH THIS PROGRAM IS USUALLY RELATED TO IMPROPER 
GEOMETRY DEFINITION. 


The pre -exponential factor, N= 1, has units of (cm-particle-sec). The 
activation energy, N= 3, has the units (cal/mol). 
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SAMPLE CASE FOR FINITE RATE CHEMISTRY 
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PROGRAM SEAGULL 

SUPERSONIC INVISCID FLOW ANALYSIS 
NASA LANGLEY RESEARCH CENTER 
RELEASED 02/01/76 
LAST MODIFICATION 02/25/77 


TWO DIMENSIONAL FLOW 

TWO DIMENSIONAL TEST CASE WITH SHOCK + INTERACT 1 
OUTPUT AT EVERY 1 STEPS 

GAMMA - • 1 230+C1 DISTANCE: .9000+01 TOLERANCE: . 1000+02 SCALE 

ICN-~ 1 IBN: 1 ISO: G ISY: 0 I T Y : 0 ITP: C I VI: 0 

.2000+02 DEGREES BREAK IN THE UPPER WALL AT Z: .1500+01 

.1000+02 DEGREES PREAK IN THE LOWER WALL AT Z- .1000+01 


INITIAL CONDITIONS 


TEMPERATURE (OEG 

. KELVIN) 

1 .0000000+03 

VELOCITY (M/S) 


5.3167825+03 

PRESSURE(ATM) 


6.8050000-01 

OENSITT (G/CC) 


1 .1281265-04 

MOLECULAR WT. 


1 . 3602172+01 

gamma 


1.2846758+00 

TOTAL ENTHALPY 


5.7722905+11 

ENTROPY 


2.1301588+09 

MOLE FRACTION 

H 

3.4969999-02 

MOLE FRACTION 

H 2 

2.3927000-01 

MOLE FRACTION 

H 20 

6.9387999-01 

MOLE FRACTION 

0 

2. 4700000-03 

MOLE FRACTION 

OH 

2.6460000-02 

MOLE FRACTION 

02 

2 . 94OOCQQ-03 

MOLE FRACTION 

H02 

0.0000000 


<D 


FINITE RATE CHEMISTRY CALCULATIONS 


4 > 

O 


H2 - 32 C OHBUST ION KI NETICS 


REACTIONS BEIN6 CONSIDERED 


1 

H 

♦ 

H 

♦ 

M2 

I 

H2 

♦ 

M2 

~Z 

H 

♦ 

0 

♦ 

Ml 

r 

OH 

♦ 

Ml 

■3 

H 

4 

OH 

♦ 

Ml 

z 

H20 

♦ 

Ml 

4 " 

H ‘ 

♦ 

02 

♦ 

M3 

z 

H02 

4 

M3 

5 

0 

♦ 

0 

♦ 

Ml 

z 

02 

4 

Ml 

--- 

OH 

4 

H 



z 

H2 

4 

0 

7 

OH 

4 

0 



z 

H 

4 

02 

- 8 

OH 

4 

H2 



I 

H20 

4 

H 

9 

OH 

♦ 

OH 



z 

H20 

4 

0 

itr 

HO 2 

♦ 

H 



z 

H 2 

4 

02 

u 

H02 

♦ 

H 



z 

OH 

4 

OH 

12 

H02 

4 

H2 



z 

H20 

4 

OH 

13 

H02 

♦ 

0 



z 

OH 

4 

02 

14 

H02 

♦ 

OH 



z 

H20 

4 

02 


KFr**EXP(B/RT**M)/T**N 




A 

N 

B 

H 

2 

2 

1 . 016+1 8 

1.00' 

.0 

.00 

2 

1 

3.628+1 5 

. DO 

.0 

.00 

2 

2 

2. 213+22 

2.00 

.0 

.00 

2 

3 

1.995+15 

.00 

1000.0 

.00 

2 

4 

1. 379+18 

1.00 

-340.0 

.00 

1 

4 

8.432+09 

-1.00 

-7000.0 

.00 

1 

1 

2. 4D9+1 3 

.00 

. 0 

.00 

1 

4 

6. 023+06 

-2. '00 

-2900.0 

.00 

1 

3 

6.023+12 

.00 

-1100.0 

.00 

1 

3 

2.530+13 

.00 

-700.0 

.00 

1 

3 

2.530+14 

.00 

-1900. 0 

.00 

1 

3 

7.228+11 

.00 

-18700. 0 

.00 

1 

3 

4. 819+13 

.00 

-1000. 0 

.00 

1 

3 

4.999+13 

.00 

-1000.0 

.00 


CATALYTIC SPECIES BEING CONSIDERED 


Ml 

= 1.00 

H 

, 1.00 

H2 

, 6.00 

H20 

, 1.00 

0 

, 1.00 

OH 

, 1.00 

02 

r 

1.00 

H02 

V 

M2 

-20.00 

H 

, 2. so 

H2 

.10.00 

H20 

, 1.00 

0 

, 1.00 

OH 

, 1.00 

02 

4 

1.00 

H02 

t 

M3 

= 1.00 

H 

, 2.50 

H2 

,16.00 

H20 

t 1.00 

0 

, 1.00 

OH 

, 1.00 

02 

f 

1.00 

H02 

t 


K= I Z- .9000+00 D2= .0000 


COJTETT tfALC -FORCE'S: 

FR=" .0000 

FZ = 

.0000 Mr 

.0000 

VFR r 

.0000 VFZr 

.0300 

VHr .0000 


UPPER 

U ALL FORCES 

FR = • 0000 

FZ = 

.0000 M= 

.0000 

VFR r 

.0000 VFZ= 

.0000 

V Mr .0000 


Stm T)F 

FORCES 

.0000 


.0000 

.0000 


.0000 

.0000 

.0000 


FORCE 

ON COMPUTATIONAL PLANE = 

.0000 

AREAr . 

0000 






— HOKFNT 

-ABOUT R AX] 

. 0003 ~ 

HOME NI 

ABOUT Z AXISr .0000 






— r 


-R — ' 

v 

tr 

4 

T 

S 

M 

THE ' 

RHO 

HTOT 

c 

21 

.1000+01 

. 1000+01 

.0000 

.6801+01 

. 1000+01 

.0000 

.6000+01 

.0000 

. 1000+01 

.1538+01 

-c- 

"20“ 

— . 94 7 4+trtr - - 

_ .inoo+m — 

.0000 

.6801+01 

. 1000+01 

.0000 

.6000+01 

.0000 

.1000+01 

.1538+01 

0 

19 

.8947+00 

. 1000 + 01 

.0000 

.6RU1+01 

.1000*01 

.0000 

.6000+01 

.0000 

. 1000+01 

.1538+01 

.... 0. 

- rs- 

2T423+TI0— 

-.1000+01 

.0000 

.6801+01 

.1000+01 

.0000 

.6000+01 

.0003 

. 1000*01 

.1538+01 

L 

17 

.7895+00 

. 100C+01 

.0000 

.6801+01 

.1000+01 

. DODO 

.6000+01 

.0000 

. 1000+01 

.1538+01 

o 

■- rr 

. TT6 8+OD 

- .1000+01 

.0000 

.6801+01 

.1000+01 

.0000 

.6000+01 

.0000 

.1000+01 

.1538+01 

0 

15 

.6842+00 

.1000+01 

.0000 

.6801+01 

. 1000+01 

.0000 

.600 0+01 

.0000 

. 1000+01 

.1538+01 

Cr 

-1-4- 

.631 6+00 

.1000*01 

.0000 

.6801+01 

.1000+01 

. 0000 

.6000+01 

iOOOO 

.1600+01 

.1538+01 

0 

1 3 

.5789+00 

.1000+01 

.0000 

.6801+01 

.100G+01 

.0000 

.6000+01 

.0000 

. 1000+01 

.1538+01 

c 

TT 

.5263+00 " 

.1000+01 

-.0000 

.6801+01 

.1000+01 

.0000 

.6030+01 

.0000 - 

-.1000+01 - 

i 1538*01 

c 

1 1 

.4737+00 

• 1000 + 01 

.0000 

.6801+01 

.1000+U1 

.0000 

.6030+01 

.0000 

. 1000+01 

.1538+01 

o— nr 

— .-*711*00 

.1000+01 

,0000 

.6801+01 

.1000+01 

.0000 

.6030+01 

.0000 

.1000+01 

.1538+01 

0 

9 

.3684+00 

.1000+01 

• ooco 

.6801+01 

. 100Q+C1 

.0000 

.6000+01 

.0000 

. 1000+01 

.1538+01 



c 

'8 

. 3153*00 

. 1000*01 

.0000 

. 68 C 1 +0 1 

. 1000+01 

.0000 

. 6000*01 

.0000 

. 1000+01 

. 1538+01 

0 

7 

. 2632*00 

. 1000+01 

.0000 

. 6801+01 

. 1000*01 

.0000 

. 6000+01 

.0000 

. 1000+01 

. 1538+01 

0 

6 

. 2106*00 

. 1000*01 

.0000 

. 6801+01 

. 1000+01 

.0000 

. 6000+01 

.ooon 

. 1000+01 

. 1538+01 

G 

5 

. 1579*00 

. 1000+01 

.ocoo 

. 6801+01 

. 1000+01 

.0000 

. 6000 + C 1 

.0300 

. 1000+01 

. 1538+01 

0 

4 

. 1053*00 

. 1000+01 

.0000 

. 6801+01 

. 1000*01 

.0000 

. 6000+01 

.0000 

. 1000+01 

. 1538+01 

0 

3 

. 5263-01 

. 1000+01 

.0000 

. 6801+01 

. 1000*01 

.0000 

. 6000+01 

.0000 

. 1000+01 

. 1538+01 

“ O ' 

- _ 2 . 

.0000 " 

. 1000*01 

.0000 

. 6801+01 

. 1000+01 

.0000 

. 6000+01 

.0000 

. 1000+01 

. 1538+01 


PT 

R ( N ) 

H 

H2 

HOLE FRACTIONS 
H20 0 

21 

1.00000 

3 . 49704-02 

2 . 39272-01 

6 . 93887-01 

2 . 47002-03 

20 

.94737 

3 . 49704-02 

2 . 39272-01 

6 . 93987-01 

2 . 47002-03 

19 

.89474 

3 . 49704-02 

2 . 39272-01 

6 . 93887-01 

2 . 47002-03 

— T8 — 

-■; 84?n 

— 7 '. 4 * 704 -02 

2 • 37272-01 

6 . 93887-01 

2 . 47002-03 

17 

.78947 

3 . 49704-02 

2 . 39272-01 

6 . 93887-01 

2 . 47002-03 

— T6 - 

. 73684 ” 

3 . 49704-02 

2 . 39272-01 

6 . 93887-01 

2 . 47002-03 

15 

.68421 

3 . 49704-02 

2 . 39272-01 

6 . 9 3 8 B 7-01 

2 . 47002-03 

14 

.63158 

3 . 49704-02 

2 . 39272-01 

6 . 93887-01 

2 . 47002-03 

13 

.57895 

3 . 49704-02 

2 . 39272-01 

6 . 93887-01 

2 . 47002-03 

- ~ T2 — 

• S263 ? 

“ 7 . 4 * 704-02 

7 .- 37272-01 

6 ”. 95787=01 

2 .' 47002-03 

11 

.47368 

3 . 49704-02 

2 . 39272-01 

6 . 93887-01 

2 . 47002-03 

10 “ 

' V4 2105 

3 . 49704-02 

2 . 39272-01 

6 . 93887-01 

2 . 47002-03 

9 

.36842 

3 . 49704-02 

2 . 39272-01 

6 . 93887-01 

2 . 47002-03 

8 ‘ 

' .31579 

3 . 49704-02 

2 . 39272-01 

6 . 93887-01 

2 . 47002-03 

7 

.26316 

3 . 49704-02 

2 . 39272-01 

6 . 93887-01 

2 . 47002-03 

‘ -6 — 

- V210S3 

3 . 49704-02 

2 . 39272-01 

6 . 93887-01 

2 . 47002-03 

5 

.15789 

3 . 49704-02 

2 . 39272-01 

6 . 93887-01 

2 . 47002-03 

4 

.10526 

3 . 49704-02 

2 . 39272-01 

6 . 93887 - C1 

2 . 47002-03 

3 

.05263 

3 . 49704-02 

2 . 39272-01 

6 . 93887-01 

2 . 47002-03 

2 

.00000 

1 . 49704-02 

2 . 39272-01 

6 . 93887-01 

2 . 47002-03 


OH 

02 

H02 

HOL UT 

GAHMA 

2 . 64603-02 

2 . 94003-03 

0.00000 

9 . 99990-01 

1 . 28467*00 

2 . 64603-02 

2 . 94G03-03 

0.00000 

9 . 99990-01 

1 . 28467+00 

2 . 64603-02 

2 . 94003-03 

0.00000 

9 . 99990-01 

1 . 26467+00 

2 . 64607-02 

2 , 94003-03 

0. OOOOO 

9 . 99990-01 

1 . 28467+00 

2 . 64603-02 

2 . 94003-03 

0.00000 

9 . 99990-01 

1 . 28467+00 

2 . 64603-02 

2 . 94003-03 

0.00000 

9 . 99990-01 

1 . 28467+00 

2 . 64603-02 

2 . 94003-03 

0.00000 

9 . 99990-01 

1 .28467 + 00 

2 . 64603-02 

2 . 94U03 - Q3 

0.00000 

9 . 99990-01 

1 . 28467+00 

2 . 64603-02 

2 . 94003-03 

0.00000 

9 . 99990-01 

1 . 28467+00 

■ 2 . 64607-02 

2 . 94003-03 

0.00000 

9 . 99990-01 

1 . 26467+00 

2 . 64603-02 

2 . 94003-03 

0.00000 

9 . 99990-01 

1 . 28467+00 

2 . 64603-02 

2 . 94003-03 

0 .00000 

9 . 99990-01 

1 . 28467*00 

2 . 64603-02 

2 . 94003-03 

0.00000 

9 , 99990-01 

1 . 28467*00 

2 . 64603-02 

2 . 94003-03 

0.00000 

9 . 99990-01 

1 . 26467+00 

2 . 64603*02 

2 . 94003-03 

0.00000 

9 . 99990-01 

1 . 28467+00 

2 . 64603-02 

2 . 94003-03 

OiOOOOO 

9 . 99990-01 

1 . 26467+00 

2 . 64603-02 

2 . 94003-03 

0.00000 

9 . 99990-01 

1 . 28467*00 

2 . 64607-02 

2 . 94003-03 

0.00000 

9 . 99990-01 

1 . 28467+00 

2 . 64603-02 

2 . 94003-03 

0.00000 

9 . 99990-01 

1 . 26467+00 

2 . 64603-02 

2 , 94003-03 

0.00000 

9 . 99990-01 

1 . 28467+00 


•t- 

H 


SAMPLE CASE FOR VIBRATIONAL NON- EQUILIBRIUM 


-O 

S3 

iELT.OlL DATA 

ELT007 RLIB70 12/01-19:12:01-1*0) 


000001 

000 

8X0T SEAGULL 








000002 

000 

_L 100 1 l — 

0 

X? 






000003 

000 

12 

24 4 2 1> 

0 







000004 

000 

CCcC jO c a 

S E 

fc Hr* 





000005 

000 

4 • o 

1.2 430 

'.O 


1C.C 

3 . 0 




000006 

000 

434. 









000007 

000 

434 * 









000008 

000 

1 .o 









000009 

000 

-30.0 









000010 

000 

C.*98 









000011 

000 

3wj. 









000012 

000 

2.0 

2000. 4340. 







000013 

000 

0.00322 

0.00162 0.00876 

0 

.01083 0 

.01018 

0,00873 0 

.00393 


000014 

000 

0.37433 

0.07024 0.05909 

0 

.02363 0 

.42544 




000015 

000 

C02100 

44.01 -90.081 








000016 

000 

0. 

6.955 

37. 

932 

-2.074 

50. 

6.955 

37.932 

-1.725 

000017 

000 

100. 

6.955 

42. 

753 

-1.378 

150. 

6.955 

45.574 

-1.030 

ooooie 

000 

200. 

6.955 

47. 

575 

-0.682 

250. 

6.955 

49.127 

-0.334 

000019 

000 

300. 

6.955 

50. 

395 

0.013 

400. 

6.955 

52.396 

0.709 

000020 

GOO 

500. 

6.955 

53. 

948 

1.404 

600. 

6.955 

55.216 

2.100 

000021 

000 

700. 

6*955 

56. 

288 

2.796 

800. 

6.955 

57.217 

3.491 

000022 

000 

900. 

6.955 

58. 

036 

4.187 

1000. 

6.955 

58.769 

4.882 

000023 

000 

1200. 

6.955 

60. 

037 

6.273 

1400. 

6.955 

61.109 

7.664 

000024 

000 

1600. 

6.955 

62. 

038 

9.055 

1800. 

6.955 

62.857 

10.446 

000025 

000 

2000. 

6.955 

63. 

590 

11.838 

2200. 

6.955 

64.253 

13.229 

000026 

000 

2400. 

6.955 

64 . 

858 

14.620 

2600. 

6.955 

65.415 

16.011 

000027 

000 

2800. 

6.955 

65. 

930 

17.402 

3000. 

6.955 

66.410 

18.793 

000028 

000 

C020Q1 

44.01 -87.331 







000029 

000 

0. 

6.955 

37. 

932 

-2.074 

50. 

6.955 

37.932 

-1.725 

000030 

000 

100. 

6.955 

42. 

753 

-1.378 

150. 

6.955 

45.574 

-1.030 

000031 

000 

200. 

6.955 

47. 

575 

-0.682 

250. 

6.955 

49.127 

-0.334 

000032 

000 

300. 

6.955 

50. 

395 

0.013 

400. 

6.955 

52.396 

0.709 

000033 

000 

500. 

6.955 

53. 

948 

1.404 

600. 

6.955 

55.216 

2.100 

000034 

000 

700. 

6.955 

56. 

288 

2.796 

800. 

6.955 

57.217 

5.491 

000035 

000 

900. 

6.955 

58. 

036 

4.187 

1000. 

6.955 

58.769 

4.882 

000036 

000 

1200. 

6.955 

60. 

037 

6.273 

1400. 

6.955 

61.109 

7.664 

000037 

000 

1600. 

6.955 

62. 

038 

9.055 

1800. 

6.955 

62.857 

10.446 

000038 

000 

2000. 

6.955 

63. 

590 

11.836 

2200. 

6 . 955 

64.253 

13.229 

000039 

000 

2400. 

6.955 

64 . 

.858 

14.620 

2600. 

6.955 

65.415 

16.011 

000040 

000 

2800. 

6.955 

65. 

930 

17.402 

3000. 

6.955 

66.410 

18.793 

000041 

000 

C 02000 

44.01 -94.051 









000042 

000 

Q. 

6.955 

37.932 

-2.1)74 

50. 

6.955 

37.932 

-1.725 

000043 

oco 

ion. 

6.955 

42.753 

-1.378 

150. 

6.955 

45.574 

-1.030 

000044 

000 

200. 

6.955 

47.575 

—0.68? 

250. 

6.955 

49.127 

-0.334 

000045 

000 

300. 

6.955 

50.395 

n . a i 3 

400. 

6.955 

52.396 

0.709 

000046 

coo 

500. 

6.955 

53.943 

1 .404 

60n. 

6.955 

55.216 

2.100 

000047 

coo 

700. 

6.955 

56.286 

2.796 

80C. 

6.955 

57.217 

3.491 

000045 

000 

900. 

6.955 

58.036 

4.187 

1000. 

6.956 

S8.769 

4 .88? 

000049 

liOO 

1200. 

6.955 

60.037 

6.273 

1400. 

6.955 

61.109 

7.664 

000050 

000 

1600. 

6.955 

62.038 

9.055 

1800. 

6.955 

62.857 

10.446 

000051 

□ 00 

2000. 

6.955 

63.590 

11.838 

2200. 

6.955 

64.253 

13.229 

000052 

000 

2400. 

6.955 

64.856 

14.620 

2600. 

6.955 

65.415 

16.011 

000053 

000 

2800. 

6.955 

65.930 

17.402 

3000. 

6.955 

66.410 

18.793 

000054 

000 

C02110 44.01 

-88. 

111 






000055 

000 

0. 

6.955 

39.31 

-2.074 

50. 

6.955 

39.31 

-1.725 

000056 

uoo 

100. 

6.955 

44.131 

-1.378 

150. 

6.955 

46.951 

-1.030 

000057 

000 

200. 

b. 955 

48.952 

-0.68? 

250. 

6.955 

50.504 

-0.334 

000058 

000 

300. 

6.955 

51 .77Z 

0.013 

400. 

6.955 

53.773 

0.709 

000059 

000 

500. 

6.955 

55.325 

1.404 

600. 

6.955 

56.593 

2.100 

000060 

000 

7uQ ■ 

6.955 

57.665 

2.796 

BOO. 

^>.955 

58.594 

3.491 

000061 

uoo 

900. 

6.955 

59.413 

4.187 

1000. 

6.955 

60.146 

4.882 

000062 

000 

1200. 

6.955 

61.414 

6.273 

1400. 

6.955 

62.487 

7.664 

000063 

000 

1600. 

6.955 

63.415 

9.055 

1800. 

6.955 

64.235 

10.446 

000064 

000 

2000. 

6.955 

64.967 

11.836 

2200. 

6.955 

6 5 • 6 3 

13.229 

000065 

coo 

2400. 

6.955 

b6 • 2 35 

14.620 

2600. 

6.955 

66.792 

16.011 

000066 

000 

2800. 

6.955 

67.308 

17.402 

3000. 

6.955 

67.788 

18.793 

000067 

uoo 

C02010 44.01 

-92. 

141 






000068 

uoo 

0. 

6.955 

39.31 

-2.074 

5n. 

6.955 

39.31 

-1.725 

000069 

□ 00 

100. 

6.955 

44.131 

-1.376 

150. 

6.955 

46.951 

-1.030 

00007G 

000 

2 00. 

6.955 

48.952 

-0.682 

250. 

6.955 

50.504 

-0.334 

000071 

UOO 

300. 

6.955 

51.772 

0.013 

400. 

6.955 

53.773 

0.709 

000072 

000 

500. 

6.955 

55.325 

1.404 

600. 

6.955 

56.593 

2.100 

U00073 

000 

700. 

6.955 

57.665 

2.796 

800. 

6.955 

58.594 

3.491 

000074 

000 

900. 

6.955 

59.413 

4.187 

1000. 

6.955 

60.146 

4.882 

000075 

000 

1200. 

6.955 

61.414 

6.273 

1400. 

6.955 

62.487 

7.664 

000076 

uoo 

1600. 

6.955 

63.415 

9.055 

1800. 

6.955 

64.235 

10.446 

000077 

000 

2000. 

6.955 

64 .967 

11.838 

2200. 

6.955 

65.63 

13.229 

000078 

ooc 

2400. 

6.955 

66.235 

14.620 

2600. 

6.955 

66.792 

16.011 

000079 

000 

2800. C 

6.955 

67.306 

17.402 

3000.0 

6.955 

67.788 

18.793 

000060 

□ 00 

C02Q20 44.01 

-90. 

284 






000051 

uoo 

0. 

fa. 955 

40.115 

-2.074 

50. 

6.955 

40.115 

-1.725 

000062 

000 

100. 

6.955 

44.936 

-1.378 

150. 

6.955 

47.757 

-1.030 

ooooej 

ooc 

200. 

6.955 

49.758 

-0.682 

250. 

6.955 

51.310 

-0.334 

000064 

oco 

300. 

6.955 

52.578 

0.013 

400. 

6,955 

54.579 

0.709 

000065 

□ GO 

500. 

6.955 

56.131 

1 .4 04 

600. 

6.955 

57.399 

2.100 

000066 

COO 

700. 

6.955 

58.471 

2.796 

8 00. 

6.955 

59.400 

3.491 

000067 

oou 

900. 

6.955 

60.219 

4.187 

1000. 

6.955 

60.952 

4.88? 

000068 

COO 

1200. 

6.955 

62.220 

6.273 

1400. 

6.955 

63.292 

7.664 

000089 

OQU 

1600 . 

6.955 

64.221 

9.055 

1800. 

6.955 

64.040 

10.446 

CQU090 

000 

2000. 

6.955 

65.773 

11.838 

2200. 

6.955 

66.436 

13.229 

000091 

GOO 

2400. 

b.955 

66.04 1 

14.620 

2600 . 

6.9 55 

67.598 

16.011 

00009? 

uoo 

2800. 

6 • V 5 5 

68.113 

17.402 

3000. 

6.955 

63.593 

18.793 

000093 

000 

C0203Q 44.01 

-88. 

531 







OOOOS 4 

GOO 

0. 


6.955 

QQ0Q95 

000 

100. 


6.955 

000096 

GOO 

200. 


6.955 

000097 

GOO 

300 . 


6.955 

000098 

000 

5 GO . 


6.955 

000099 

GOO 

700. 


6.955 

Q0Q1U0 

000 

900. 


6.955 

00010! 

JCO 

1200. 


6.955 

000102 

00 u 

1 6 uO . 


6.955 

000103 

000 

2000. 


6.955 

000104 

oou 

2400. 


6.955 

□00105 

000 

281)0. 


6.955 

000106 

000 

N?(0> 

28.016 

0. 

000107 

000 

0. 


6.955 

000108 

uoo 

100. 


6.956 

000109 

000 

200. 


6.957 

000110 

000 

300. 


6.958 

000111 

000 

5 bo * 


6.960 

000112 

GOD 

700. 


6.965 

000113 

000 

900. 


6.972 

0001 14 

uoo 

1200. 


6.983 

000115 

000 

1600. 


6.999 

000116 

oou 

2000. 


7.015 

000117 

GOO 

2400. 


7.G31 

000118 

000 

2800. 


7.045 

000119 

GOO 

N 2 { 1 ) 

26.016 

6, 

000120 

000 

0. 


6.955 

000121 

GOO 

100. 


6.956 

000122 

GOO 

200. 


6.957 

000123 

GOG 

300. 


6.958 

000124 

QOO 

500. 


6.960 

000125 

000 

700. 


6.965 

000126 

000 

900. 


6.972 

000127 

000 

1200. 


6.983 

000128 

OOU 

1600. 


6.999 

000129 

000 

2 000. 


7.015 

000130 

OOU 

2400. 


7.031 

000131 

000 

2800. 


7.045 

000132 

ODD 

02 

31.9988 0 

000133 

000 

0. 


6.953 

000134 

000 

loo. 


6.958 

000135 

000 

200. 


6.961 

000136 

000 

3 L)Q * 


7.023 

000137 

000 

500. 


7.431 

00013fl 

000 

700. 


7. 863 

0001 39 

000 

9 UO . 


8.212 

U0D14D 

000 

1200. 


8.52 7 

000141 

coo 

1600 . 


0.900 

000142 

uoo 

2300. 


9.029 

000143 

oou 

2400. 


9.248 

000144 

noo 

2800. 


9.465 

00D145 

000 

H20 18.016 

>57. 


40.687 

-2.074 

50. 

45.506 

-1 . 37B 

150. 

50.329 

-0.662 

250. 

53.150 

0.013 

400. 

56.703 

1.404 

600. 

59.043 

2.796 

BOO , 

60.791 

4.187 

1000. 

62.792 

6.273 

1400. 

64.793 

9.055 

1800. 

66.345 

11.838 

2200. 

67.613 

14.620 

2600. 

68.665 

17.4U2 

3000. 

3.796 

71,4 


33.339 

-2.075 

5 0. 

38.160 

-1.379 

150. 

42.932 

-0.663 

250. 

45.803 

0.013 

400. 

49.357 

1.405 

600. 

51.700 

2.797 

800. 

53.451 

4.190 

1000. 

55.456 

6.284 

1400. 

57.470 

8 . 9 b 3 

1800. 

59.033 

11.883 

2200. 

60.313 

14.692 

2600. 

61.399 

17.507 

3000. 

3.798 

71.4 


33.339 

-2.075 

50. 

38 .160 

-1.379 

ISO. 

42 .982 

-0.663 

2 50. 

45.803 

0.013 

400. 

49.357 

1.405 

600. 

51.700 

2.797 

800. 

53.451 

4.190 

1000. 

55.458 

6.284 

1400. 

57 .470 

8.963 

1800. 

69.033 

11.863 

2200. 

60.313 

14.692 

2600 . 

61.399 

17.507 

3000. 

36.572 

-2.075 

50. 

41.395 

-1.381 

150. 

46 .218 

-0.635 

250. 

49.047 

0.013 

400. 

52.7 22 

1.455 

600. 

5S.297 

2.988 

600. 

57.320 

4.600 

1000. 

59.729 

7.114 

1400. 

62.222 

10.563 

1800. 

64 .210 

14.149 

2200 • 

65.876 

17.804 

2600. 

67.317 

21.545 

30U0 


6.955 

40.687 

-1 .725 

6.955 

48.329 

-1.030 

6.955 

51.882 

-0.334 

6.955 

55.151 

0.709 

6.955 

5 7. 971 

2.100 

6.955 

59.972 

3.491 

6.955 

61.524 

4 .382 

6.955 

63.664 

7 .664 

6.955 

65.612 

10.446 

6.955 

67.008 

13.229 

6.955 

68.170 

16.011 

6.955 

69.165 

18.793 

6.955 

33.339 

-1.727 

6.956 

40.981 

-1.031 

6.957 

44.535 

-D.335 

6.959 

47.805 

0.709 

6.962 

50.627 

2.101 

6.968 

52.631 

3.494 

6.975 

54.186 

4.888 

6.991 

56.535 

7.681 

7.008 

58.294 

10.481 

7.02? 

59.70? 

13.28 7 

7.038 

60.876 

16.099 

7.052 

61.885 

18.917 

6.955 

33.339 

-1.727 

6.956 

40.981 

-1.031 

6.957 

44.535 

-0.335 

6.959 

47.805 

0.709 

6.962 

50.627 

2.101 

6.968 

52.631 

3.494 

6.975 

54.186 

4.888 

6.991 

56.535 

7.681 

7.008 

53.294 

10.481 

7.022 

59.702 

13,26 7 

7.038 

60.876 

16.099 

7.052 

61.885 

18.917 

6.958 

36.572 

-1.729 

6.960 

44.216 

-1.033 

6.992 

47.769 

-0.337 

7.196 

51.091 

0.724 

7.670 

54.090 

2,210 

8.063 

56. 3bl 

3*786 

8.336 

58.192 

5.427 

0.674 

61.055 

8.835 

8.916 

63.265 

12.354 

9.139 

65.076 

15.966 

9.354 

66.620 

19.664 

9.551 

67.973 

23.446 


000146 DOG 0.0 □ . 0 0.0 -2. 36? 50. 3.9805 18.198 -1.974 

000147 000 100. 7.961 36.396 -1.561 150. 7.965 39.156 -1.16 

D00148 000 200. 7.969 41.916 - .764 250. 7.993 43.535 - .384 

000149 „ COO 300. 6.027 45.155 .015 400. 8.186 47.433 .825 

000150 000 500. e.415 49.334 1.654 600. S.676 50.891 2.509 

000151 000 700. 8.954 52.249 3.390 800. 9.246 53.464 4.300 

000152 00 U 900. D 9.547 54.570 5.240 1000.0 9.851 55.591 6.209 

000153 000 1200.0 10.444 57.440 8.240 1400.0 10.987 59.092 10.384 

000154 000 1600.0 11.462 60.591 12.630 1800.0 11.869 61.965 14.964 

000155 000 2000.0 12.214 63.234 17.373 2200.0 12.505 64.412 19.846 

000156 000 2400.0 12.753 65.511 22.372 2600.0 12.965 66.540 24.945 

000157 000 2800.0 13. 14b 67.508 27.556 3000.0 13.304 68.420 30.201 

000158 COO HE 4.003 0. 


000159 

000 

0. 


4.968 21.255 

-1.481 

50. 


4.968 

21.255 

-1.233 

000160 

000 

100. 


4.968 24.698 

-0.934 

150. 


4 .968 

26.128 

-0.736 

000161 

000 

200. 


4.968 28.142 

-0.488 

250. 


4.968 

29.048 

-0.239 

000162 

000 

300. 


4 , 9o8 30.156 

0.009 

400. 


4.968 

31.586 

0.506 

Q001b3 

COO 

500. 


4.968 32.694 

1.003 

600. 


4.968 

33.600 

1.500 

000164 

000 

700. 


4.968 34.366 

1.996 

BUO . 


4 ,966 

35.029 

2.493 

0DO165 

000 

900. 


4.968 35.614 

2.990 

1000. 


4.968 

36.138 

3.487 

000166 

000 

12 00. 


4 • 9t>8 37.044 

4.481 

14 00., 


4.968 

37.809 

5.4 74 

000167 

000 

1600. 


4.968 38.473 

6.468 

1800. 


4.968 

39.058 

7.461 

000168 

coo 

2000. 


4.968 39.581 

8.455 

2200. 


4.968 

40. OSS 

9.449 

000169 

000 

2400. 


4.968 40.437 

10.442 

2600. 


4.968 

40.885 

11.436 

000170 

000 

2800. 


4.968 41.253 

12.430 

3000. 


4.968 

41.596 

13.423 

0001 71 

000 

hi 









□001 72 

000 

1.0 1.0 

1.0 

1.0 1.0 1.0 1.0 

0.0 0 

.0 0.0 

i • 

0 0 . n 1 

.0 


000173 

000 

M2 









000174 

000 

1.0 1.0 

1.0 

1.0 l.Q 1.0 1.0 

1.0 1 

.0 1.0 

1. 

0 1.5 1 

.0 


000175 

000 

M3 









000176 

000 

1.0 l.U 

1.0 

1.0 1.0 1.0 1.0 

0.5 0 

.5 0.6 

1. 

0 0.0 1 

.0 


0001 77 

000 

M4 









000176 

000 

0.0 0.0 

0.0 

0.0 O.U 0.0 U.O 

2.0 2.3 2.4 

0. 

0 1.0 0 

.0 


□00179 

000 

C02 

GQL 

REACTION MOUEL (AIR 

-BREATHING) 19 

OCT 1976 26 

REACTIONS 


□00180 

000 

C021 10+C02000 

=C02100+C02010 



12 

1 .25-13-0 

.5 


000181 

000 

C02Q30+C02000 

-C02 100+C02Q 10 



12 

1.8 -15-0 

.5 


000182 

□ 00 

C02030+C02000 

-C02020*C02010 



12 

3.1 -13-0 

.5 


000133 

DOO 

C02100+CG2000 

;C02010*C02010 



11 

4.0 -13 



000164 

000 

C02020+C02000 

-C0201Q+C02010 



12 

1.4 -12-0 

.5 


000185 

000 

C020U1 *K1 


ICO 2 1 10«-Hl 



64 

1.1 -27-4 

•8 1484. 


000186 

000 

C02GC1 *M4 


-C02 1 1 0+M4 



64 

1.9 -31-5 

.8 2436. 


0001S7 

uno 

COZOUl ♦Hi 


=C02U30*M1 



64 

8.1 -31-5 

•6 1484. 


0001 38 

□ CO 

C02001 +M4 


=C02030+M4 



64 

1.4 -34-6 

.6 2436. 


000189 

000 

C021 10+M2 


=C02030+M2 



62 

4.3 -17-1 

.5 


000190 

oou 

C021 10+K2 


-C02020+M2 



64 

4.5 -27-4 

.2 - 903. 



-P- 

Ln 


-IN 

ON 


000191 

000 

C02110+M2 

-002020*^2 

64 

8.8 -20-2.5 

- 4410. 

000192 

uoo 

C021 10+M2 

-CQ? 1 00*M2 

64 

8.6 -24-3.8 

- 549. 

000193 

000 

C02030+M2 

=C02020+M2 

64 

9.3 -22-3.3 

- 1230. 

000194 

000 

C02030+M2 

=C02100*M2 

64 

1.1 -21-3.0 

- 1060. 

000195 

000 

C02100+M2 

-C02G20+M2 

62 

8.0 -18-1.5 


000196 

000 

C02100+M3 

- C020 1 0* H 3 

64 

5.6 -22-3.3 

- 1480. 

000197 

000 

C021 OO+HE 

-C02010+HE 

14 

1.8 -21-3.0 

843. 

000196 

000 

C02020 ♦ M3 

-C02010+M3 

64 

2.1 -21-3.2 

- 1350. 

000199 

000 

C02020+HE 

=C0201G*HE 

14 

3.8 -21-3.0 

843. 

000200 

000 

C0201 0*M3 

-C02000+H3 

64 

3.4 -26-4.2 

1130. 

000201 

000 

C020 1 0 *HE 

-C02000*HE 

14 

9.9 -22-3.0 

843. 

000202 

000 

C020Q1+N2 (0) 

=C02000+N2( 1 > 

12 

8.3 -12 0.5 


000203 

000 

C02010+H20 

-C02000+H20 

15 

3.2 -13 

45 

000204 

uoo 

C02001 +H20 

=C02030+H20 

11 

4.0 -13 


000205 

uoo 

N2(l) ♦H20 

-N2 ( 0 ) + H20 

15 

1.1 -10 

139 

ENQ ELT. 








SADD f P DATA 


0.333 
0. 333 


dXQT SEAGULL 



PROGRAM SEAGULL 

SUPERSONIC INVISCID FLOW ANALYSIS 
NASA LANGLEY RESEARCH CENTER 
RELEASED Q2/C1/76 
LAST MODIFICATION 02/25/77 


TWO DIMENSIONAL FLOW 
CCCC 0000 CASE *** 
OUTPUT AT EVERY 100 STEPS 
GAMMA: . 1200*01 DISTANCE: 

ICiN: 1 ISN: 1 ISO: u ISY: 
-.3000+02 OEGREES BREAK IN 
-.3000+02 DEGREES BREAK IN 


.4300+03 TOLERANCE: 

0 I T Y : 0 ITP: U IVI 
THE UPPER WALL AT Z: 
THE LOWER WALL AT Z: 


.1000+02 SCALE: .3000+01 

0 

. 1000+01 

.9980+00 


INITIAL CONDITIONS 


TEMPERATUREtDEG. 

KELVIN) 

2 • OQOuOOQ +P 3 

VELOCITY (M/S) 


4.5628436+03 

PRESSURE ( ATM ) 


2.0000000+00 

DENSITY (G/CC ) 


2.2614246-04 

MOLECULAR WT. 


1.8554989+01 

GAMMA 



1.4520176+00 

TOTAL 

ENTHALPY 


2. 1616025+12 

ENTROPY 


2. 1521133+09 

MOLE 

FRACTION 

C021CD 

3.2200000-03 

MOLE 

FRACTION 

C02GD1 

1.6200000-03 

MOLE 

fraction 

C02Q00 

8.7599999-03 

MOLE 

FRACTION 

C0211D 

1 . 083U000-02 

MOLE 

FRACTION 

C0201 0 

1.0180000-02 

MOLE 

fraction 

C02D20 

8.7299999-03 

MOLE 

fraction 

C 02 0 3 0 

3.9299999-03 

MOLE 

FRACTION 

N2(0> 

3.7432999-01 

MOLL 

FRACTION 

N2( 1) 

7.0239999-02 

MOLE 

fraction 

02 

5.9090000-02 

MOLE 

fraction 

H 2 0 

2.3630000-02 

MOLE 

FRACTION 

HE 

4.2543999-01 


■p- 

00 


C02 o 


OL REA 


CHON 


REACTIONS BEING CONSIDERED 


1 

C021 10+ 

C 0 2000 

z 

CO 2 1 00+ 

C02010 

2 

C02030 ♦ 

C02000 

: 

C02100+ 

C02010 

3 

C02030+ 

CO 20 0 0 

: 

C02D20+ 

C02010 

4 

C02100+ 

C02000 

= 

C02Q10+ 

C02013 

5 

• C02D20+ 

C02Q0G 

: 

CO 20 1 0 ♦ 

C02010 

6 

C02001+ 

Ml 

r 

C02110+ 

Ml 

7 

C02001+ 

M4 

= 

C021 10+ 

H4 

a 

C02001+ 

Ml 


C02030+ 

Ml 

9 

C02001+ 

M4 

: 

C02030+ 

M 4 

10 

C021 10+ 

M2 

: 

C02030+ 

M2 

1 1 

C 02 1 1 0 ♦ 

M2 

: 

C02020+ 

M2 

12 

C021 10+ 

M2 

z 

C02020+ 

M2 

13 

C021 10+ 

M2 

z 

C02100+ 

M2 

14 

C02030+ 

H2 

z 

C 02020 + 

M2 

15 

C02030+ 

M2 

z 

C02100+ 

M2 

16 

C021G0+ 

M2 

z 

CO 2020 ♦ 

M2 

17 

C02100+ 

M3 

z 

C02010+ 

M3 

18 

C02100+ 

HE 

z 

CO 20 1 0 + 

HE 

19 

C02020+ 

M3 

z 

C02010+ 

M3 

20 

C02020+ 

HE 

z 

C02010+ 

HE 

21 

C02010+ 

M3 

z 

CO 2000 * 

M3 

22 

C02010+ 

HE 

z 

C02000+ 

HE 

23 

C02001 ♦ 

N2 (0) 

z 

CO 2000+ 

N2 ( 1 ) 

24 

C02010+ 

H20 

z 

CO 2000 ♦ 

H20 

25 

C02Q01 ♦ 

H20 

: 

C02030+ 

H 20 

26 

N2(l) ♦ 

H20 

: 

N2 ( 0 ) ♦ 

H20 


CATALYTIC SPECIES BEING CONSIDERED 



Ml 

= 

1 .00 

CO 2 1 0 Ot 

1 

.00 

C020U1 , 

1.00 

C02000 




.CO 

N2 (0 1 , 


.00 

N2m , 

.00 

02 


M2 

= 

1 . 00 

C021UU, 

1 

.00 

C0230 1 • 

1 .GO 

CO 2000 




1.00 

N2<01 , 

1 

.00 

N2(l) , 

1.00 

02 


M3 

= 

1.00 

C02100. 

1 

.00 

C0200 1 $ 

1.00 

C02000 




.50 

N2(0) , 


.50 

N2 ( 1 > , 

.60 

02 


H4 

= 

.00 

CO 2 100 t 


.00 

CO 200 1 , 

.00 

CO 2000 




2.00 

N2 ( 0 ) , 

2 

.00 

N2 ( 1 ) , 

2.40 

02 

K: D 2: 

.9950+00 

02 = 

.0000 






LOWER WALL 

FORCES 

FR: 

.OUOO F2 

r 


.0000 

M= 

.0000 

upper WALL 

FORCES 

FR = 

.0000 F 2 

: 


.0000 

M = 

.0000 

SUM OF FORCES 


.OCOO 



.0000 


• 0000 


FORCE ON COMPUTATIONAL PLANE: .0000 AREA: .OOGC 

MOMENT ABOUT R AXIS: .0000 HOHENT ABOUT Z AXIS: .00 


MODEL 


< AIR-B 


REATHI 


NG| 1 


KF:A*LXPtB/RT**M»/T**N 
AN B H 


1 2 7.529+10 -.50 

1 2 1.084+09 -.50 

1 2 1.867*11 -.50 

1 1 2.409*11 .00 

1 2 8.432*11 -.50 

6 4 6.626-04 -4.80 

6 4 1.144-07 -5.80 

6 4 4.879-07 -5.60 

6 4 8.432-11 -6.60 

6 2 2.590*07 -1.50 

6 4 2.710-03 -4.20 

6 4 5.300*04 -2.50 

6 4 5.180*00 -3.80 

6 4 5 .602*02 -3.30 

6 4 6.626*02 -3.00 

6 2 4.819*06 -1.50 

6 4 3.373*02 -3.30 

1 4 1,084*03 -3.00 

6 4 1.265*03 -3.20 

1 4 2.289+03 -3.00 

6 4 2.048-02 -4.20 

1 4 5.963*02 -3.00 

1 2 4.999+12 .50 

1 5 1.927*11 .00 

1 1 2.409+11 .00 

1 5 6.626+13 .00 


1.30 C02 1 1 0 t 1*00 C02010 
1.00 H20 , .00 HE 

1.00 CO 2 1 1 0 1 1.00 C02010 
1.00 H20 , 1.50 HE 

1.00 C02 1 1 0 t 1.00 C02010 
1.00 H20 , .00 HE 

.00 C02110, .00 C02Q10 

.00 H20 , 1.00 HE 


VFR: .0000 VF2= 

VFR: .0000 VFZ: 

.0000 


.0 .00 

.0 .00 

.0 .00 

.0 .00 

.0 .00 

1484.0 .00 

2436.0 .00 

1484.0 .00 

2436.0 .00 

.0 .00 

-903.0 .00 

-4410.0 .00 

-549.0 .00 

-1230.0 .DO 

-1060.0 .00 
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